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Abstract 


MINOS  is  a Fortran  system  for  solving  large-scale  linearly 
constrained  optimization  problems.  The  System  Manual  gives  an 
overview  of  the  system,  the  programming  conventions  used,  data  struc- 

I 

I tures,  tolerances,  and  error  conditions.  Details  are  given  of  a 

I practical  implementation  of  the  method  of  Bartels  and  Golub  for 

! maintaining  a sparse  LU  factorization.  The  reduced-gradient  approach 

for  handling  a nonlinear  objective  function  has  been  described 
elsewhere  by  Murtagh  and  Saunders;  further  implementation  details 
are  included  here.  The  System  Manual  should  facilitate  interfacing 

of  MINOS  with  other  optimization  software. 
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1.  INTRODUCTIOM 


MINOS  is  an  in-core,  Fortran-based  optimization  system  for  the 
solution  of  large-scale  linear  and  nonlinear  programming  problems  involving 
sparse  linear  constraints.  Full  user  information  is  given  in  the  MINOS 
User's  Guide  (Murtagh  and  Saunders,  1977).  Theoretical  aspects  of  the 
algorithm,  and  experience  with  an  earlier  version  of  the  system,  are  des- 
cribed in  Murtagh  and  Saunders  (1978). 

This  manual  provides  details  of  the  implementation  methods  used.  Much 
of  the  complexity  of  the  system  is  concerned  with  a reliable  and  efficient 
implementation  of  the  primal  simplex  method,  as  described  in  Saunders  (1976). 
The  nonlinear  progranroing  features  are  an  integral  part  of  the  design,  and 
this  required  that  most  of  the  standard  simplex  routines  be  moderately 
generalized.  However,  the  bulk  of  the  code  for  nonlinear  problems  is  con- 
tained in  24  additional  subroutines. 

Occasional  reference  is  made  to  MINOS/GRG,  an  extension  of  the  present 
system  designed  to  handle  problems  with  nonlinear  constraints  and  a sparse 
Jacobian  matrix  (Jain,  Lasdon  and  Saunders,  1976).  This  system  will  be 
documented  elsewhere. 

The  primary  goals  throughout  the  coding  of  MINOS  have  been  as  follows 
(in  alphabetical  order). 

Efficiency  - Use  of  efficient  algorithms;  conservation  of  core  require- 
ments through  dynamic  storage  allocation;  avoidance  of 
multi-dimensional  arrays;  sequential  access  to  data  where 
possible  to  avoid  frequent  page-swapping  in  a virtual 
memory  environment. 

Flexibility  - Free-format  input  for  user-defined  parameters;  ability  to 
solve  problems  of  arbitrary  size  without  recompilation  of 
the  source  program;  various  restart  facilities;  provision 
for  interface  with  existing  commercial  systems. 

Modularity  - Well-defined  functions  for  each  subroutine;  use  of  para- 
meters to  pass  all  array  storage  and  dimension  information. 
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Portabi 1 i ty 


Readability  - 


Use  of  only  those  Fortran  constructs  that  are  commonly 
available  on  (large)  machines;  separation  and  internal 
documentation  of  routines  that  are  necessarily  machine- 
dependent;  calculation  of  precision  tolerances  from  one 
machine-dependent  constant. 

Indentation  of  loops  as  a bare  minimum;  meaningful  names 
for  variables;  a reasonable  amount  of  in-line  documenta- 


Reliability  - Use  of  numerically  stable  algorithms;  consistency  checks 
on  input  data,  restarting  files  and  gradient  computation. 

There  is  always  room  for  improvement  in  putting  such  goals  into 
practice,  but  certainly  some  such  aims  are  necessary  if  a system  is  to 
develop  while  the  algorithms  being  implemented  are  themselves  continually 
evolving. 


I 

2.  OVERVIEW 

2.1  Storage  Allocation 

To  ensure  storage  economy,  a single  array  Z(*)  is  used  to  provide  all 
array  workspace.  The  type  of  Z must  match  the  largest  word  size  used  (e.g. 

REAL*8  for  IBM,  REAL  for  CDC).  Arrays  of  other  types  are  allocated  compactly 
within  Z,  using  three  constants  to  define  the  word  sizes  available  on  any 
particular  machine. 

Array  storage  is  made  available  to  all  subroutines  through  parameters 
rather  than  through  blank  or  labelled  COMMON.  (Use  of  COMMON  diminishes  the  j 

effectiveness  of  optimizing  compilers.)  It  should  be  noted  that  sensible 
modularization  eliminates  potential  inefficiency  in  this  respect.  If  a sub-  j 

routine  performs  any  significant  work  at  all,  or  of  it  is  called  infrequently,  I 

the  overhead  in  passing  20  or  30  parameters  is  negligibly  more  than  in  passing 
just  one. 

A particular  advantage  arising  from  the  use  of  a single  array  Z is  the 
ability  to  acquire  core  according  to  problem  size  without  recompiling  the 
program.  For  the  IBM  version  of  MINOS,  an  assembly  language  main  program  has 
been  provided  by  Hedges  (1975),  allowing  core  requirements  to  be  defined  at 
run-time  through  the  JCL  REGION  parameter.  An  analogous  Algol  procedure  is 
available  for  the  Burroughs  B6700  version,  which  incidentally  circumvents  the 
restriction  in  B6700  Fortran  that  arrays  be  no  longer  than  65535  words.  This 
is  documented  in  Saunders  (1977). 


2.2  Subroutine  Hierarchy 

The  diagram  below  illustrates  the  broad  grouping  of  subroutines.  Once 
a problem  has  been  input,  the  dimensions  of  most  arrays  are  determined. 

Hence  subroutine  MINOS  is  able  to  pass  many  segments  of  the  array  Z to 
DRIVER,  which  thereafter  knows  them  as  separate  arrays  with  meaningful  names. 
A few  other  arrays  for  the  LU  file  remain  to  be  allocated  by  INVERT  within 
the  remaining  part  of  Z. 


Input 


Solve  and  Output 


A more  complete  classification  of  subroutines  is  given  on  the  next  page. 
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Input  routines  (for  processing  the  SPfCS  and  MPS  files) 


HASH  HOVE  BPS  HPSIh  SPECS  SPECS2 


Basis  loading  routines 


CPASH  IHSEHT  LOADB  LOADN  NHSPCH 


Basis  and  solution  saving  routines 


ECJMP  POMCH  SAVE?  SCIN  SOLPRT  STATE 


Algorithwic 

routines 

required 

for  purely 

linear  prograes 

BTBANl 

TOPHC 

ITEPOP 

PRICE 

BTPASO 

FTPAHL 

LPITH 

P3 

EHKFS 

FTFANU 

HINOS 

P4 

CHUZF 

GO 

HKLIST 

SETPI 

DEIVEF 

INITLZ 

HODLU 

SETX 

PACTOF 

INVEFT 

PACKLU 

TPNSVL 

OHFACK 


Additional 

algorithiic 

routines 

for  nonlinear  probleis 

ADDCCL 

CHKGPD 

FrjNGRD 

B1  ADD 

ALIGN 

CHOZC 

GETGPD 

BIHOD 

CAICFG 

CCMDFP 

NEBPIC 

B 1EF0D 

CAICG 

PEICOL 

BESEIH 

E1SUE 

CG 

DOT 

RGITN 

SEARCH 

CHKDIP 

DOT! 

RTFSOL 

SQUEEZ 
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3.  PROGRAMMING  CONVENTIONS 


The  originating  source  code  is  for  the  IBM  Systems  360  and  370.  The 
notes  here  will  facilitate  conversion  to  other  machines. 

* Types 

In  most  subroutines  the  type  of  a variable  is  defined  by  the  first 
character  of  its  name,  as  follows: 

! A-B  REAL 

t 

C-G,  0-Z  REAL*8 

H INTEGER*2 

I-N  INTEGER 

The  only  variables  beginning  with  H are  half-word  arrays.  These  and  other 
arrays  are  always  typed  and  dimensioned  explicitly,  e.g. 

INTEGER*2  HA(NE) 

REAL  A(NE) 

REAL*8  X(MN) 

' A few  variables  are  explicitly  typed  LOGICAL. 

1 IMPLICIT 

Non-standard  implicit  typing  is  used  only  for  double  precision  simple 
l variables,  and  is  accomplished  where  necessary  by  the  statement 

[ IMPLICIT  REAL*8(C-G,0-Z) 

I 

[ (which  immediately  follows  a SUBROUTINE  statement). 

COMMON 

Blank  COMMON  is  not  used.  There  are  24  labelled  COMMON  blocks,  con- 
taining only  simple  variables  and  a few  very  short  arrays.  None  of  the 
COMMON  variables  is  used  for  variable  dimensions. 

' EQUIVALENCE 

No  EQUIVALENCE  statements  are  used,  except  in  a straightforward  manner 
in  subroutine  MPS: 
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EQUIVALENCE  (KEY,  ID(1),  ID2(1),  ID3(1),  ID5(1)) 


DATA 

Sone  subroutines  use  data  statements  to  initialize  integer  variables 
to  strings  of  at  most  four  characters,  e.g. 

DATA  LRHS  /3HRHS/ 

Such  items  immediately  precede  the  first  executable  statement  in  a sub- 
routine. For  CDC  compatibility,  apostrophes  are  never  used. 


Variable  dimensions 

Subroutine  parameters  are  used  to  define  the  size  of  all  variable-length 
arrays.  Arrays  are  strictly  one-dimensional.  Declarations  such  as  REAL  A(l) 
are  never  used. 


Subscripts 

Subscript  expressions  of  the  form  X(M+J)  are  often  used,  but  more 
complex  forms  such  as  X(HB(J))  are  not. 


Avoidance  of  trouble-spots 

The  following  Fortran  syntactical  items  are  not  used: 

arithmetic  IF 
assigned  GO  TO 
BACKSPACE 
blank  COMMON 
BLOCK  DATA 
COMPLEX 
ENDFILE 
ENTRY 

EQUIVALENCE  (except  as  noted  above) 

EXTERNAL 

FIND 

FUNCTION 

NAMELIST 

nonstandard  RETURN 
PAUSE 
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PRINT 

PUNCH 

TIME 

unformatted  READ  or  WRITE 
$ in  identifiers 
/ around  subroutine  parameters 

3. 1 Conversion  to  Other  Machines 

(aj  INTEGER*2  may  be  replaced  by  INTEGERbb  throughout  (where  b means 
"space"),  although  most  compilers  do  the  right  thing  anyway. 

(b)  If  double  precision  is  not  required,  the  following  modifications 
should  be  made: 


change 

to 

in 

REAL*8 

REALbb 

all 

DABS 

ABS 

all 

DMAXl 

AMAXl 

subroutine  INITLZ 

DMINl 

AMINl 

subroutine  SEARCH 

DSQRT 

SORT 

all 

D- 

E- 

subroutine  INITLZ 

DLOG 

ALOG 

subroutine  CALCFG 

The  IMPLICIT  statements  may  then  be  deleted. 

(c)  On  the  Burroughs  B6700,  a simpler  alternative  is  to  use  the  compiler 
control  card  $SET  DBLTOSNGL. 

(d)  To  speed  editing  it  may  help  to  note  that  the  specifications 
IMPLICIT 
INTEGER*? 

REAL*8 

all  begin  in  column  7. 

(e)  The  quantities  EPS,  NWORDR,  NWORDI,  NWORDH  must  be  set  correctly  in 
subroutine  INITLZ. 

(f)  Floating-point  constants  are  coded  in  single  precision,  e.g. 

PLINFY  = l.OE+30,  and  therefore  do  not  require  conversion.  (An 
exception  is  the  arguments  to  DMAXl  in  subroutine  INITLZ,  requiring 
the  change  of  D-  to  E-  as  noted  above.) 

10 
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[ (g)  If  double  precision  is  required  but  the  IMPLICIT  statement  is  not 

I allowed,  all  relevant  simple  variables  beginning  with  C-G,  0-Z  must 

I be  typed  explicitly,  e.g.  DOUBLE  PRECISION  EPS,EPS0,EPS1.  This 

\ includes  local  variables,  COMMON  variables  and  subroutine  parameters. 

[ Such  a conversion  will  require  many  additional  lines  of  otherwise 

^ unnecessary  code.  It  will  also  be  highly  prone  to  error,  since 

[ Fortran  compilers  normally  do  not  issue  warnings  or  error  messages 

^ for  undeclared  variables. 

‘ (h)  For  CDC-RUN,  several  READ  statements  of  the  form 

READ(u,f,END=m)  1 

must  have  the  ",END=m"  eliminated,  and  additional  cards  inserted  of 
the  form 
IF  (E0F,u)  m,n 
n CONTINUE 

(i)  In  FORMAT  statements,  apostrophes  are  used  to  delimit  strings.  How- 
ever, strings  never  contain  the  character  *.  Hence  for  CDC  machines 
it  is  acceptable  to  convert  every  apostrophe  to  * (except  perhaps  in 
comment  cards). 
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4.  VARIABLES 

4.1  Varidtle-diaension  irfcrnaticr, 

The  follcwlng  IHTE6EE  vaciatlds  aie  used  to  define  the  current  length 
of  various  var iable**dirensicn  arrays.  They  usually  occur  at  the  front 
of  a subroutine's  paraieter  list. 

LFN  Length  of  hash  table  used  to  store  tow  names  during  input. 

K No.  cf  rows  in  constraint  matrix  A. 

MAXL  "laxiaum  no.  cf  elements  currently  allcwed  for  in  L file. 

KAXS  Maximum  no.  of  superhasic  variables  allowed. 

MAXU  Maximum  no.  of  elements  currently  allowed  for  in  U file. 

MAXZ  Size  of  Z,  the  single  array  of  cere  in  which  all  else  resides. 

MCOIS  No.  of  columns  allowed  for  in  A during  input  (must  be  >=  N\  . 

HELMS  No.  cf  nonzeros  allcwed  for  in  A during  input  (must  be  >=  'i?)  . 
MN  d ♦ HAXS.  HaximuH  dimension  of  basics  superhasics. 

KNN  max(  HN,  NN  ) . Used  only  for  divension  of  work  vectors  Y and 

n. 

MPONS  No.  of  rows  allcwed  for  in  A during  input  (must  be  >=  Hi . 

HS  H ♦ NS.  Current  dimension  of  basics  ♦ superbasics. 

HI  M ♦ 1. 

H2  H ♦ KIN?.  Up  to  this  many  tra.nsfotmations  allowed  for  in  f.  and 

0. 

N No.  of  variables  In  constraint  matrix  (structurals  ♦ PHS  ♦ 

slacks)  . 

NE  No.  of  elements  in  constraint  matrix. 

SEN  No.  of  entries  in  hash  table  during  input. 

NF  Size  of  array  F = MSPK*  (BSFKt 1) /2  where  HSPK  = max  no.  of 

spikes. 

NN  No.  of  nonlinear  variables  (those  involved  in  nonlinear 

objective) . 

NNO  ■ax(  NN,  1 ) 

NP1  Nel. 

NF  Size  of  array  F = M AXS*  (HAX.St  1) /2. 

NS  Current  no.  of  superbasic  variables. 

12 


1 


iS 


I 

i 

1 


i 


4.2  COMMON  d*»clarations 


Th«  followirg  labelled  COttNCN  blocks  are  us#>d  to  segment  paraBiatars 
into  various  logical  sets.  This  is  to  avoid  repetition  of  or.o  large 
CCNMON  block  in  every  sutrcutice. 

I' 

I 

! 

LOGICAL  CCNV,FI.gTHT  i 

COMMON  /BGCCM  / TCL SWP , NE UM P , NSPlKf , HSF K | 

COMMON  /CGCCM  / C GBET  A , ITMCG , .1 3GCG , SOPCG,  BESTPT  j.j 

COMMON  /CONVCM/  I T A £ H ,ET  A FG  ,C  C N V ( 4)  , L VLTCL  , NXTPHE  , NF  A 1 1 [; 

COM.MON  /C03E  / KZ1,K72,KZJ  i 

COMMON  /DJCna  / TCLDJ  1,TOLDJ2,TCLDJ3,TOLDJ  ji 

COMMON  /EPSCCM/  E E S , EESC , EPS  1 , E t S2 , EPS  3 , EPS4 , EPS5 , PLIN F Y | 

COMMON  /FILES  / I SC F , IN F UT , 10 LC E , I NEW P , I NSET , IPNC H , ILO AD , I DUMP  1 

COMMON  /FFEQS  / KC  H F , KI K V , K S A V , KIOG  , 1 1 F PEC  , 1 2FBE0 

COMMON  /FXCOM  / E X , SI  NF  , KTCE  J , MI  NT  B Z , Nf  X , N IN  F , ICB  J , NPE  CIJ  I: 

COMMON  /INTCCH/  IT N , ITN  LI F , NP H S , KMCE LU  , K MODP I j; 

COMMON  /INVCCM/  I N VFC » IN V IT N , IN VMOL , NBELEM , K PL, KDL , KHU , KL U , KFF  [ 

COM.MON  /ITNLCG/  D J Q,  T HE T A , PI  VOT  ,C0  K D , NT  NOPT  , JP , JQ , K APD , KSU B [ 

CCMHON  /ITNLG2/  I H E AD  , .IC  1 , 0Q2  , ■1':  1 , J F 2 

COMMON  /LPCO.M  / KR  H S , NS  1 , « A XF  , IF  FP  , I DEP  UG  i 

COMMON  /LUFILE/  1 B 3F C , S FT A L , N 2? A P , NET AU , KL , K U, LPGN , JLEGN 
COMMON  /HPSCCM/  AI J TOI , E ST E UC  ( 2)  , M I N MA X , MLST , ME F , 

1 NAME  (2)  ,MOEJ  (2)  ,MRHS  (2)  ,MFNG  (2)  ,rf3ND  (2) 

COMMON  /NLCCCM/  M N L FC N , L t F El V , LPRI KI 
COMMON  /PAPMC.M/  D P A PM  ( 9)  , T P AP  B (O) 

COMMON  /PPCCC  M/  N E A F P F , NPRC , JPPC, K P PC , M FE J ,N FE J , J RE J ( 20) 

COMMON  /PRCCM2/  N B U IF F , NEfc SB 

COMMON  /RGTCLS/  X"-C  L ( 3)  , F'^C  I ( 3)  , GTC  I ( 3)  , PINOFM  , PG NOFM  , TO!.  E G 
COMMON  /SCLNCM/  I F Cl N , K ' C I N , M SC  I N , NST AT F , LCHKGB 
COMMON  /TOLS  / TC L X , TC I PI  V , TE PI V 1 , TB P J V 2 , TOLPO W , X NOFM 
COMMON  /WOEDSZ/  N fcOF DF , NWOI DI , NWCBDH 
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4.3  reflnltion  of  CCH.^CN  variables 


AIJTCL 
BSTFUC  (2) 

CGBETA 

COMD 

CCNV  (4) 

DJO 

DPAPH  (9) 

EPS 

BPSO 

EPS1 
BPS  2 
EPS  3 
EPS4 
EPS5 
*TARG 

ETASH 
PTOL  (3) 

?X 

GTCL  (3) 

IDEBOG 

ICUMF 

TEPP 

IBEAE 

IlOAD 

IMBED 


Tolerance  for  ignoring  matrix  coefficients  during  MPS  input. 

Default  lower  and  upper  bound  values  for  structural 
variables. 

Conjugate  gradient  scalar,  (new  D)  = -gradient  ♦ CGDETA*(oli 
D)  . 

Square  of  condition  number  of  diagonal  of  P.  This  provides 
a lower  bound  on  the  cendition  of  the  Hessian  approximation. 

Logical  variables  showing  result  cf  convergence  tests  in 
PGITH. 

Best  reduced  cost  found  by  PPICE. 

Array  of  floating-point  parameters. 

Machine  precision.  1fi**(-13)  for  IBM  370. 

EPS**(4/5).  Used,  among  ether  things,  for  packing  L'l 
tranforaatiens. 

EPS**  (2/3).  Used  fer  finite  differencing  in  ADDCOL. 

SCPT  (EPS)  . 

EPS**( 1/3) . 

EPS**(1/4)  . 

EPS**(1/5) . 

Hinimization-within-sntspaca  accuracy  tolerance,  used  in 
PGITN. 

Linessarch  accuracy  tclerar.ee,  used  in  MEWPTC. 

Current,  loose  and  tight  tcls  or.  change  in  objective. 

Current  objective  value  (linear  ♦ nonlinear  terns) . 

Current,  loose  and  tight  tols  on  nern  cf  reduced  qradi* 

Debug  level. 

Unit  no.  for  DUMP  file. 

Error  indicator  for  numerous  routines. 

Determines  if  ancther  beading  is  to  be  printed  in  itr.  log. 
Unit  no.  for  LOAD  file. 

Determines  if  part  of  L should  be  imbedded  ir  A. 
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IMENE 

INPUT 

INSat 

IKTITN 

INVMOO 

INVFt 

lOBJ 

lOLCE 

IPASK  (9) 

lENCH 

ISCP 

ISOLN 

ITN 

ITNCG 

ITNLi:' 

I1FFF9 

I2PFE0 

JLBGN 

JP 

JPEC 

JC 

J01 

JQ2 

JBEJ  (20) 

JP1 

JF2 

KADD 

KCHK 


Unit  no.  for  KFW  BIT-P.AP  file. 

Unit  no.  for  EFS  file. 

Ur.it  no.  foe  INSIFT  filf-. 

No.  of  iterations  performeJ  since  last  invert. 

No.  of  modifications  sade  to  LU  factors  sir.ee  last  invert. 
Invert  requast.  Set  positive  for  nuierous  reasons. 

Row  no.  of  linear  objective.  zero  if  none. 

Unit  no.  for  CID  BIT-MAF  file. 

Array  of  inteqer  parameters. 

Unit  no.  for  F'JNCM  file. 

Unit  no.  for  SCRATCH  file. 

Unit  no.  for  SCIUTTCN  file. 

No.  of  iterations  since  start  of  run. 

Vo.  of  consecutive  cen  juqata-qradicnt  itr.s  since  restart, 
aaxinum  nc.  of  iterations  allowed. 

May  be  set  ty  I1F9E0  card  in  SPECS  file.  1 causes  the  bumo 
and  spike  pattern  to  he  displayed  after  each  INVERT. 

May  be  set  ty  I2FEE0  card  in  SPFCS  file. 

Beginning  in  HL,CL  cf  transf ermatien  pcintad  to  by  LBGN. 

The  JP-th  ccluran  of  the  basis  will  be  replaced. 

Colunn  no.  marking  end  of  last  partition  of  A scannad  by 
PRICE. 

Column  no.  of  variable  selected  by  PRICE  or  CHUZO  to  ‘'nter 
basis . 

Column  no,  of  variable  entering  ♦■he  basis  (for  printim). 

Column  no,  of  variable  entering  the  superbasic  set. 

Column  DOS.  of  variables  rejected  frem  basis  by  INVERT. 
Column  no,  of  variable  leaving  the  basis. 

Column  nc.  of  variable  leaving  the  superbasic  set. 

Indicator  for  modi f ica tio.n  of  R. 
check  frequency  (testing  tow  residuals). 


r 


KDL  Beginning  of  array  DI  in  aeoiory  array  Z.  (Part  of  LU 

KDO  Beginning  cf  array  DU  ir  Bemory  arra  7,  (Part  of  III  file.) 

KFF  Beginning  of  array  F in  Beuoty  array  7.  (Part  of  LU  fil#*.) 

KHL  Beginning  of  array  HI  in  memory  array  7i.  (Part  of  Ltl  file.) 

KHD  Beginning  of  array  HU  in  memory  array  Z.  (Part  of  LU  tile.) 

KINV  Factorization  (invert)  freguency. 

KL  First  free  position  in  I file,  kl-1  = no.  of  nonzeros  in  L. 

KLOG  leg  frequency  (iteration  log). 

KMODLU  Set  to  1 if  lu  factors  are  to  be  modified, 

KMODII  set  to  1 if  FI  is  to  be  reccpputed. 

KPBC  Hhich  partition  cf  A mas  last  scanned  ty  PPICF. 

KPHS  Column  no.  of  the  PHS  (treated  as  a variable  fixed  at  -1.0). 

KSAV  Save  frequency  (for  PEH  BIT  PAP). 

KSOLH  Indicator  for  final  call  to  CAICFG  at  optimal  solution. 

KSUB  Indicator  fet  modification  cf  P. 

KU  First  free  position  in  0 file.  KO-1  = no.  of  nonzercs  in  U. 

This  includes  nonzeros  in  columns  of  U that  havo  been 
deleted  during  updating. 

KZ1  Points  to  beginning  cf  core  available  for  INVEST. 

KZ2  Not  used. 

KZ3  Points  to  end  of  core  used  by  arrays  other  than  for  LU. 

LPGN  The  first  L transformation  to  be  used  by  FTPANL  and  BlfANL. 

LCHKGR  Indicator  for  verifying  gradient  of  oblective. 

LDEPIV  Derivative  level.  Not  used. 

LPFINT  Print  level.  Net  used. 

LVLTOL  2 (3)  means  loose  (tight)  tcls  for  normal  (final) 

iterations . 

MAXH  Maximum  dixensicn  of  P.  Max  no.  of  superbasics  for 

quasi-Nemton. 

nPND(2)  Name  of  BCUND  set. 

MEF  RaxifflUB  no.  of  error  messages  (cf  each  type)  during  MPS 

input . 

MINIMZ  1 for  min,  -1  for  max. 
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HINNAX 
MLST 
«CBJ  »2) 
HODCG 
HFEJ 

BPHS  (2) 
KFEG  (2) 
MSGCG 


f»SCLN 

MSPK 

NAME  (2) 

NSELEM 

NEUME 

NETAL 

NPTAF 

NETAU 

NEWSE 

NFAII 

NFX 

NINF 

NMaiEN 

NNLFCN 

NONCFT 

NPARF? 

NPHS 

N?EC 

NPFCE 

NPFJ 


Contains  relevant  string  'MIN'  cr  'MAX'. 

laxiBUB  lines  cf  MFS  liata  to  ba  lipte'l  during  input. 

Naae  of  OBJECTIVE  tow. 

Defines  which  ccr.  jugate-gradient  alaorithm  is  to  bo  used. 

Maiiaura  no,  of  rejected  vars.  to  be  saved  for  PANDAID  in 
PPICE. 

Naae  of  flHS. 

Name  of  BRNGE  set. 

Used  for  detecting  first  entry  to  conjugate-gradient  routine 
(CG)  . 

Indicator  of  request  to  print  solution. 

Maximum  nc.  of  spikes  currently  allowed  for. 

Naae  of  problem  (or  first  card  cf  MPS  file). 

No.  cf  nonzetoE  in  basis  cn  entry  to  INVEST. 

No.  of  bumps  feund  by  P4  during  refactorization. 

Ho.  of  transformations  in  L after  invert. 

No.  of  transf craatiors  added  to  L during  updating. 

No.  of  transformations  in  U,  includina  ones  flagged  as 
deleted . 

No.  of  nontasic  variables  selected  to  becoae  superbasic. 

Nc.  of  consecutive  linesearch  failures. 

No.  of  evaluatiens  of  objective  function  and  its  gradient. 
No.  of  infeasibilities. 

Multiple  price  indicator.  Overrides  NFA^PP. 

No.  of  nonlinear  functions  (0  or  1 in  MINOs.  More  ir.  '''PG)  . 

No.  of  nonoptimal  columns  during  last  PPICEing  sweep. 

Partial  price  indicator  = no.  of  partitions  into  which  the 
censtraint  matrix  A is  divided  for  the  purposes  of  oricing. 

Current  phase  (1,  2,  3 cr  U) . 

Size  of  partitions  of  A for  partial  pricing. 

Problem  no.,  for  use  in  CALCFG. 

No.  cf  variables  rejected  from  basis  by  INVERT. 


■j 

J 
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NSPIFE 


Ho.  of  spikes  in  F (active  spikes  in  D) 
P4,F»CTCP,BCDITJ. 


set  by 


i 


NSTATE 

NS1 

NHOFCH 

NHOBDI 

NHOBDB 

NXTPHS 

PIHOBH 

PIVCT 

PLIHFY 

PFSTFX 

BGNOPH 

SINP 

THETA 

TOLDJ 

TC1DJ1 

TOLrj2 

TCLDJ3 

TOLPIV 

TOLRG 

TCLBCH 

TCLSHP 

TCLX 

TPPIV1 

TBP1V2 

HTOBJ 

XNOBH 
XTOL  (3) 


Paraaeter  of  CALCFG.  0,1,2  indicate  roraal,  first,  last 
entries. 

No.  of  nonlinear  variables  currently  in  the  basis. 

Ho.  of  HALF  INTEGEP  words  par  word  of  aeaory  array  Z. 

Ho.  of  IHTBGEF  wcrds  pet  word  of  meiory  array  Z. 

Ho.  of  PEAL  words  per  word  of  aieBoty  array  Z. 

In  PGITN,  suggests  HPHS  for  next  itt.  FOHNC  aay  override. 
Horn  of  PI. 

Pivot  elenant  for  coluan  entering  basis. 

1.0E*30.  Used  for  "infinite"  bounds. 

Bequests  restart  cf  conjugate  gradient  algorithm. 

Hora  of  reduced  gradient  (largest  amongst  superbasics) . 

Sub  of  infeasibilities. 

Step  aoved  in  current  search  direction. 

Current  reduced  cost  tol  (scaled  by  PINOFH) . 

Reduced  cost  tol  while  infeasible. 

Initial  reduced  cost  tcl  during  feasible  itns  if  partial 
pricing. 

Final  reduced  cost  tcl  during  feasible  itns. 

MiniBum  size  for  PIvCT  daring  iterations. 

Level  to  which  FGSCPH  lust  be  reduced  before  adding 
superbasics. 

Tol  used  In  seasutlng  "BOW  EPPOB"  (size  of  EHS  - A*X) . 

LU  MOD  TCIEFAKCE.  Eelative  pivot  tol  for  updating  LH 
factors. 


) 


Feasibility  tclerance  for 

LO  ROW  TCIFBAHCE.  Relati 

LO  COL  TOLERANCE.  Pelati 

Height  on  the  true  linear 
the  coaposite  objective  b 

Nora  of  X for  basic  varia 

Current,  loose  and  tight 
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4.4  Arrays 


All  arrays  are  aeqients  ct  the  array  Z.  Most  are  allocatad  after 
input,  once  the  dlnensions  of  the  prcblen  are  known.  Those  cot.tair.inq 
the  factoritaticn  cf  the  lasts  are  reallocated  at  each  call  to 
INVEPT. 


Array 

Length 

A 

HE 

BL 

N 

DU 

N 

C 

NNO 

DL 

MAXL 

T)U 

SAXU 

? 

NF 

G 

.HAXS 

GN’EW 

HAXS 

GPD 

HN 

G1 

NHO 

HA 

HE 

HP 

MN 

HS 

NP1 

HL 

HA  XL 

HPIVF 

P2 

HPIVl 

"2 

liFIVO 

.*12 

HS  N 


Nontarc  coefficients  in  the  constraint  matrix 
( \ BHS  I ),  listed  cclumn-wise. 

Lower  bounds. 

Oppet  fccunds. 

Gradient  cf  ncrlinear  part  of  objective  function. 
Nonzecos  in  L t ransf oriatiors. 

Monzercs  in  '}  tran.sformations . 

A full  trianqular  ratrix,  part  of  L'J  factors  of  basis. 

Reduced  gradient  vector  for  supertastes. 

Sew  reduced  gradient  vector  for  superbasics. 

Gradient  vector  for  basics  and  supertasics. 

Teapctacy  gradient  vector  for  nonlinear  cbjective. 

Bow  indices  for  each  cclumn  cf  A. 

Column  ncs.  of  basics  and  superbasics. 

For  basics  (i  <-  d <=  f*)  , liS  (J)  is  r he  column  that 
pivots  cr  rew  j. 

HE(J)  points  to  the  end  of  column  d-1  in  arrays  A, HA. 

Bow  Indices  for  renzetes  In  each  transformation  in  L. 

The  first  index  for  each  tran  sf  cmat  Lon  is  neoat  iv^  if 
the  pivot  elesent  is  l.d,  otherwise  the  pivot  is  in  PT.. 

HPIVF  (d)  is  the  row  on  which  the  J-th  spike  pivots. 

HPIVL  (d)  is  the  pivot  row  for  the  J-th  transfer  mat ron 
in  I . 

HPIVU  (J)  is  the  LOW  on  which  the  J-th  column  of  il  pivots 
If  negative,  this  column  of  J has  been  deleted. 

HS (J)  is  the  state  of  the  J-th  variable  in  ( A FHH  I ) 

0 =>  variable  J is  ncntasic  at  its  lower  bound. 

1 =>  variable  j is  nonbasic  at  it.s  upper  bound. 


2 *>  Taclakle  J is  supecbasic. 

3 3>  variable  J is  tasic. 

Nontasic  IRTE  variables  may  have  state  0 or  1 


HSPIKZ 

H 

During  INVEET  defines  the  spike  structure  of  the 

basis.  (This  structure  eay  be  displayed  by  setting 
I1FBE0  =1).  During  iterations,  HSPIKE  is  used  as  a 
work  vector  in  calls  tc  ETFANl. 

HO 

HAXU 

Row  indices  for  nenzeros  in  0. 

NL 

B2 

HEL  * ML  (J)  is  the  number  of  nonzeros  in  the  J-th 
transfer oaticn  cf  I.  If  HEL  < 0,  the  transformation 
is  actually  coluwn  "-HEL"  of  A,  with  non-unit  pivot 
eleaent  pointed  to  (indirectly)  by  an  element  of  the 
array  HL. 

HO 

H2 

HO  (J)  is  the  no.  of  nenzeros  in  the  J-th  transfermation 
of  0. 

PI 

H 

The  pricing  vector  (i.e.  the  Lagrange  multipliers  or  the 
simplex  multipliers  or  the  shadow  prices  for  the  general 
constraints) . 

B 

HP 

Upper  triangular  matrix  used  to  approximate  the  reduced 
Hessian,  thus:  B'B  = Z'GZ  approximately,  where  Z 
spans  the  null  space  c£  the  active  constraints. 

X 

HH 

Values  of  the  basic  and  superbasic  variables. 

XN 

MHO 

Values  of  the  nonlinear  variables. 

T 

HNM 

Work  vector. 

11 

HMH 

Hotk  vector. 

Z 

RAXZ 

The  single  array  of  memory  in  which  all  else  resides. 
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5.  TOLERANCES 


I 


_ * 


Most  tolerances  are  defined  in  terms  of  the  machine  precision 
EPS  = e,  which  is  the  smallest  positive  number  such  that  1,0  + EPS  is 
greater  than  1.0,  when  the  arithmetic  i'^  performed  in  the  apprc  riate 
precision  (e.g.  single  on  CDC  and  Burroughs  machines,  double  on  IBM  and 
Uni  vac). 


For  convenience,  the  labelled  common  block  EPSCOM  contains  the 
following  quantities,  which  are  initialized  in  subroutine  INITLZ  and  never 
altered: 


EPS 

= 

e 

EPSO 

= 

,4/5 

EPSl 

= 

EPS2 

= 

,1/2 

EPS3 

= 

EPS4 

,1/4 

EPS5 

= 

PLINFY 

= 

10^° 

A tolerance  TOLZ  is  used  for  packing  transformations  (PACKLU,  MODLU)  and 
for  detecting  negligible  elements  in  transformed  vectors  in  order  to  speed 
transformation  processing  (BTRANL,  BTRANU,  FTRANL,  FTRANU).  This  tolerance 
is  set  to  EPSO  at  the  beginning  of  the  subroutines  mentioned.  A larger 
value  for  TOLZ  would  reduce  the  accuracy  with  which  the  constraints  Ax  = b 
could  be  satisfied,  without  ensuring  any  significant  improvement  in  storage 
or  run  time.  (In  fact,  the  effect  could  be  negative.) 

The  following  tolerances  are  set  in  subroutine  INITLZ  (see  the  listing 
of  that  subroutine  in  section  10),  The  values  used  should  be  satisfactory 
for  most  well -scaled  problems. 


Tolerance 

Value 

Description 

TOLX 

max(E^,  10"^) 

Primal  feasibility  tolerance.  Structural  and 

slack  variables  are  allowed  to  lie  outside  their 

bounds  by  as  much  as  TOLX. 

TOLDJl 

max(£^,  10"®) 

Relative  dual  infeasibility  tolerance  while 
infeasible.  In  PRICE,  the  actual  tolerance 

kC, 
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T0LDJ2 

T0LDJ3 


TOLPIV 


TOLROW 


used  to  screen  reduced  costs  (reduced 
gradients)  on  nonbasic  variables  is 
TOLDOUPINORM. 

1.0  Initial  and  final  relative  tolerances  on  reduced 

niax(e^.  10”®)  gradients  in  PRICE  when  feasible.  If  partial 

pricing  is  in  effect,  a relative  tolerance  TOLDJ 
begins  at  T0LDJ2  and  is  reduced  in  stages  to  the 
level  T0LDJ3,  Again,  the  actual  tolerance  used 
in  PRICE  is  T0LDJ*PIN0RM. 

If  partial  pricing  is  not  in  effect,  the  value 

used  in  PRICE  is  T0LDJ3*PIN0RM. 
u 

e The  smallest  pivot  allowed  in  CHUZR  during 

simplex  iterations  (LPITN).  For  non-simplex 
steps  (R6ITN)  the  value  used  is  T0LPIV*yN0RM, 
where  YNORM  is  the  norm  of  the  search  direction. 

max 10”^)  The  tolerance  used  in  SETX  to  determine  if  the 
constraints  Ax  = b are  sufficiently  well 
satisfied.  The  error  allowed  on  each  row  is 
measured  relative  to  the  norm  of  the  correspon- 
ding row  in  the  current  basis.  See  the  des- 
cription of  SETX  in  section  10. 


Certain  other  tolerances  are  available  to  the  user  to  alter  if  desired, 
via  the  SPECS  file.  (This  is  accomplished  in  subroutine  SPECS2.)  Default 
values  are  shown  below.  Fuller  details  are  given  in  the  MINOS  User's 
Guide. 


Tolerance  Default  Value 
AIJTOL  10”''° 


Description 

The  smallest  matrix  coefficient  accepted  as 
input. 


ETARG  0.2  Affects  the  accuracy  of  minimization  within 

each  subspace.  Refer  to  REDUCED  GRADIENT 
TOLERANCE. 


ETASH 


0.01  Affects  the  accuracy  of  the  linesearch. 


Refer  to  LINESEARCH  TOLERANCE. 
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TOLSWP 


TRPIVl 


TRPIV2 


WTOBJ 


1 


0,99  Relative  pivot  tolerance  for  Bartels-Golub 

updating.  Refer  to  LU  MOD  TOLERANCE  in  the 
User's  Guide. 

0.001  Relative  pivot  tolerance  for  accepting  the 

preassigned  pivot  row  in  FACTOR.  Refer  to 
LU  ROW  TOLERANCE. 

0.1  Relative  pivot  tolerance  for  accepting  one  of 

several  possible  spikes  as  substitutes  for  a 
column  whose  preassigned  pivot  row  is 
unsatisfactory.  Refer  to  LU  COL  TOLERANCE. 

0.0  The  weight  to  be  placed  on  the  true  objective 

function  when  infeasible.  (The  sum  of 
infeasibilities  always  carries  a weight  of  1.0.) 
Refer  to  WEIGHT  ON  OBJECTIVE. 


Post-script: 

TOLX  and  T0LDJ3  above  may  also  be  changed  via  the  SPECS  file. 
See  section  12. 
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6.  DATA  STRUCTURES 


Here  we  describe  the  arrays  and  variables  used  to  store  the  constraints, 
the  basis  factors,  and  the  reduced  Hessian  approximation. 


6.1  Constraint  Matrix  [A  b I] 


REAL 

A(NE) 

Nonzero  matrix  and  rhs  coefficients. 

INTEGER*2 

HA(NE) 

Corresponding  row  indices. 

INTEGER*2 

HE(NPl) 

HE(J)  points  to  the  end  of  column  J-1  in 
arrays  A and  HA. 

REAL 

BL(N) 

Lower  bounds  on  all  variables  (including  rhs 
and  slacks). 

REAL 

BU(N) 

Upper  bounds. 

Related  INTEGERS: 

M Number  of  rows  in  A (including  objective  rows). 

N Total  number  of  variables  ( structural s,  rhs  and  slacks). 

NE  Total  number  of  matrix  elements  (including  rhs  and  slacks). 

NPl  N+1. 

KRHS  Column  number  for  rhs. 

The  constraint  matrix  [A  b I]  is  stored  column-wise  in  packed  form, 
using  the  arrays  A,  HA  and  HE.  The  coefficients  for  column  J are  in  elements 
HE(J)+1  through  HE(J+1)  of  arrays  A and  HA.  The  row  indices  in  HA  for  a 
particular  column  are  not  required  to  be  in  ascending  order  (say);  they  are 
stored  in  the  same  order  as  they  occur  in  the  input  data  (which  in  general 
is  arbitrary). 

The  absolute  value  of  a particular  element  HA(I)  is  taken  as  the 
required  row  index.  Hence,  negative  entries  in  HA  may  be  introduced  in 
special  applications.  (For  example,  in  MINOS/GRG  they  are  used  to 
distinguish  variable  Jacobian  elements  from  normal  (constant)  coefficients.) 
Note  that  HA(I)=0  is  not  allowed. 

For  a typical  use  of  arrays  A,  HA  and  HE,  see  subroutine  UNPACK  in 
section  10. 


i 


I 


Explicit  upper  and  lower  bounds  are  stored  in  arrays  BU  and  BL  for  all 
variables.  The  value 

PLINFY  = 10^°  j 
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is  used  to  represent  "infinite"  bounds.  For  instance,  the  bounds  on  slack 
variables  to  implement  the  four  types  of  inequality  constraint  are  as 
f ol 1 ows : 

Constraint  Type  BL  BU 

L 0.0  PLINFY 

E 0.0  0.0 

G -PLINFY  0.0 

N -PLINFY  PLINFY 

RANGES  on  rows  are  implemented  by  altering  either  BL  or  BU  in  a trivial  way. 

Note:  It  is  more  common  for  the  constraint  matrix  to  be  stored  in  the  form 
[I  A b].  This  certainly  simplifies  the  input  procedures.  We  have 
chosen  the  less  conventional  form  [A  b I]  to  ensure  that  the 
column  numbers  of  the  first  NN  structural  variables  correspond  to 
the  nonlinear  variables  in  the  array  X of  the  user's  subroutine 
CALCFG.  This  also  allows  the  simple  test  IF  (J.LE.NN)  ...  to 
determine  whether  the  J-th  variable  is  nonlinear. 


6.2  LU  Factors  of  the  Basis 


Ignoring  certain  permutations,  subroutine  INVERT  produces  an  LU 
factorization  of  the  current  basis  B in  the  form 

li 


I u 

I 

I U 

II 

CO 

= L 

F 

— - 

F 

I 

where  L is  lower  triangular  and  F is  upper  triangular.  The  "bump  and  spike" 
structure  of  B is  vital  to  the  success  of  this  factorization.  B is  first 
permuted  to  block  lower  triangular  form  using  the  procedures  of  Duff  and 
Reid  (1975,  1976).  Each  block  ("bump")  is  then  processed  by  the  algorithm 
of  Hellerman  and  Rarick  (1971),  thereby  dividing  the  columns  of  B into 
two  sets,  called  triangle  columns  and  spike  columns.  Spikes  are  those 
columns  that  have  nonzeros  above  the  diagonal  of  the  permuted  B. 

The  LU  factorization  amounts  to  Gaussian  elimination  with  column  inter- 
changes. If  a column  in  a particular  bump  has  to  be  interchanged  with  some 
alternative  column,  the  only  alternatives  that  need  be  considered  are  the 
remaining  spikes  in  the  same  bump.  If  a triangle  column  is  interchanged 
with  a spike,  the  number  of  spikes  generally  increases  by  one.  If  two  spikes 
are  interchanged  the  number  of  spikes  generally  stays  the  same. 

The  dimension  of  F is  the  final  number  of  spikes. 


The  reasons  for  computing  the  factorization  in  thi;  way,  and  the 
(numerically  stable)  method  for  updating  L,  U and  F when  a column  of  B is 
replaced,  are  fully  described  in  Saunders  (1976).  The  result  is  a practical 
implementatidn  of  the  method  of  Bartels  and  Golub  (1969,  1971).  Our  aim  here 
is  to  describe  how  the  various  matrices  are  stored. 

For  an  initial  factorization  (immediately  after  INVERT),  the  columns  of 
L and  U are  stored  conventionally  in  packed  form,  working  forwards  through 
the  arrays  HL,  DL  and  HU,  DU.  The  upper  triangular  part  of  F is  stored  column- 
wise as  a dense  matrix  in  the  array  F. 

When  the  factorization  is  updated,  U and  F are  modified  explicitly  and 
continue  to  be  stored  in  the  same  way.  Updates  to  L are  accumulated  in  product 
form.  They  are  stored  in  arrays  HL,  DL,  much  like  columns  of  the  initial  L, 
but  interpretation  of  the  contents  of  these  arrays  is  unusually  intricate. 


Disclaimer: 


Since  the  density  of  F is  typically  about  10%  (though  in  some  cases  as 
high  as  50%),  the  storage  required  for  F becomes  excessive  if  the  number  of 
spikes  is  200  or  more.  (The  efficiency  of  FTRANU  and  BTRANU  is  not 
necessarily  degraded,  because  many  rows  and  columns  of  F can  be  skipped 
during  the  transformation  processing.) 

In  a future  version  of  MINOS  it  is  intended  to  store  F row-wise  via  a 
linked -list  data  structure. 

6.2.1  Storage  for  L 

REAL*8  DL(MAXL)  Packed  nonzeros. 

INTEGER*2  HL(MAXL)  Corresponding  indices. 

INTEGER  NL(M2)  Number  of  nonzeros  in  each  transformation. 

INTEGER*2  HPIVL(M2)  Pivot  rows. 

Related  INTEGERS: 

MAXL  Maximum  number  of  nonzeros  currently  allowed  for  in  the 

L file.  (Reset  during  each  call  to  INVERT.) 

M2  An  upper  bound  on  the  number  of  transformations  there 

could  be  before  the  next  call  to  INVERT.  (Fixed; 
M2=M+KINV,  where  KINV  is  the  refactorization  frequency.) 

NETAL  The  number  of  nontrivial  columns  in  L after  INVERT. 

(Those  corresponding  to  slacks  are  not  stored.) 

NETAR  Current  number  of  nontrivial  updates  to  L. 

KDL,  KHL,KL  See  also. 

The  K-th  transformation  in  the  L file  is  represented  by  NL(K),  HPIVL(K) 
and  a set  of  contiguous  entries  i,d  in  the  arrays  HL,DL.  Let  JL  be  a pointer 
such  that  JL+j  marks  the  j-th  pair  (i,d)  thus: 

Index 

JL 

JL+1 
JL+2 

JL+NEL 
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Not  all  items  need  be  present.  Where  relevant,  let  the  following  quantities 
be  defined: 

NEL  = NL(K) 

IPIV  = HPIVL(K)  = pivot  row 
KEY  = HL(JL) 

PIVOT  = DL(JL)  = pivot  element 

Initial  L 

For  K < NETAL,  the  K-th  transformation  represents  the  K-th  column  of  L, 
which  may  be  one  of  three  types  as  follows. 


Type  1:  NEL  < 0. 


Triangle  column  of  B imbedded  in  the 
constraint  matrix. 


This  column  of  L is  the  J-th  column  of  A,  where  J = -NEL.  Only  one 
entry  exists  in  HL  and  DL  (KEY  and  PIVOT).  KEY  is  not  used. 


Type  2:  NEL  > 0,  KEY  > 0. 


Triangle  column  of  B,  not  imbedded  in  the 
constraint  matrix. 


This  column  of  L is  the  column  vector 
(PIVOT,  d^ , 62,  ...,  dfj£L) 
with  coefficients  in  the  rows 
(IPIV,  i^,  i2,  ..., 

It  has  been  copied  explicitly  from  A into  the  L file. 


Type  3:  NEL  > 0,  KEY  < 0.  Bottom  half  of  a transformed  spike. 

This  column  of  L is  the  vector  in  Type  2 divided  by  the  pivot 
element,  i.e., 

(1.0,  di/PIV0T dN£L/PIV0T). 

The  divisions  are  not  performed  explicitly,  as  they  are  more  economic- 
ally accounted  for  whenever  the  transformation  is  used  in  FTRANL  and 
BTRANL.  Note  that  the  pivot  element  for  this  spike  transformation  is 
deliberately  chosen  to  be  1.0,  and  the  element  PIVOT  is  installed  as 
the  appropriate  diagonal  of  F.  (The  two  pivot  elements  could  be  PIVL 
and  PIVF  for  any  numbers  such  that  PIVL*PIVF  = PIVOT.  However,  the 
primary  aim  is  to  keep  L well -conditioned,  since  no  part  of  it  is  dis- 
carded, and  to  let  the  condition  of  the  basis  be  reflected  in  F.  Thus 
we  choose  L to  be  as  near  to  the  identity  matrix  as  possible.) 
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Type  2 transformations  would  not  be  required  if  triangle  columns  were 
always  imbedded  in  A.  However,  on  a machine  with  virtual  memory  and/or  a 
cache  memory,  transformation  processing  will  be  faster  with  Type  2 than  with 
Type  1.  Also  for  applications  involving  nonlinear  constraints,  the  columns 
of  A may  contain  variable  elements  of  the  Jacobian.  Their  present  values 
can  be  retained  with  Type  2,  but  not  with  Type  1. 

Updates  to  L 

For  NETAL  < K £ NETAL  + NETAR,  the  K-th  transformation  represents  the 
(K-NETAL)th  update  to  L.  The  items  KEY  and  PIVOT  do  not  exist. 

Referring  to  section  4 of  Saunders  (1976),  recall  that  any  one  such 
update  comes  from  a product  of  elementary  transformations  Ej  (j=l,2,. . . ,NEL 
say),  which  reduce  some  intermediate  matrix  F to  upper  triangular  form  F, 
using  Gaussian  elimination  with  row  interchanges.  Thus 


and 


F 


to  be  eliminated 


new  column 
of  F 


‘NEL- 


•EjE^F  = F. 


If  L is  the  previous  (possibly  updated)  factor,  the  modified  form  of  L is 

where  Mj  = Ej’^  but  since  Mj  and  Ej  are  trivially  related,  it  is  immaterial 
which  of  these  matrices  is  stored.  We  choose  to  store  Ej. 


In  the  method  of  Bartels  and  Golub,  each  transformation  Ej  is  a matrix 
of  the  form 


1 

column  i 

which  is  possibly  preceded  by  a permutation  matrix  that  interchanges  rows 
IPIVj  and  ij.  The  indices  ij  are  the  positions  of  nonzeros  being  eliminated 
from  the  bottom  row  of  F.  These  indices  and  the  corresponding  multipliers 
are  stored  as  the  pair  (ij.dj)  in  arrays  HL.DL  as  shown  earlier. 

Now,  the  intricate  part  of  the  implementation  lies  in  keeping  track 
of  the  values  IPIVj  without  storing  another  set  of  indices.  First  of  all, 
note  that  the  initial  row  index  IPIV^  is  the  pivot  row  corresponding  to  the 
column  that  was  deleted  from  the  basis.  This  value  is  stored  as  IPIV  = HPIVL(K). 

Secondly,  note  that  if  no  interchanges  were  required  to  maintain 
stability,  then  all  IPIVj  would  be  the  same  value  IPIV.  (The  update  would 
then  be  analytically  identical  to  the  method  of  Forrest  and  Tomlin  (1972).) 

In  practice  we  have  found  that 

(a)  The  number  of  eliminations  per  update  is  extremely  small.  It  is  bounded 
by  0 £ NEL  < current  no.  of  spikes,  but  typically  0 £ NEL  £ 20  because  of 
the  sparsity  of  the  bottom  row  of  F. 

(b)  It  is  frequently  the  case  that  no  interchanges  are  required,  even  if  the 
tolerance  TOLSWP  is  1.0.  (This  means  that  the  diagonal  elements  of  F 
are  usually  larger  than  the  elements  being  eliminated  from  F. ) Perhaps 
this  explains  empirically  the  practical  success  of  the  method  of  Forrest 
and  Tomlin  on  problems  that  are  not  severely  ill-conditioned. 

Because  of  point  (b)  we  distinguish  between  the  two  updating  cases. 

Type  1.  NEL  > 0.  No  interchanges. 

The  row  indices  IPIVj  are  constant,  namely  IPIV.  The  indices  ij 
are  stored  in  HL(JL+j).  Transformation  processing  in  FTRANL  and 
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I 

I 

jl 

BTRANL  can  be  coded  somewhat  more  efficiently.  In  particular,  a 
complete  transformation  can  be  skipped  in  BTRANL  if  the  element 
Y(IPIV)  in  the  vector  being  transformed  is  sufficiently  close  to 
zero. 

Type  2.  NEL  < 0.  Some  interchanges  occurred. 

An  interchange  associated  with  Ej  is  recorded  by  storing  HL(JL+j) 
negatively.  During  a forward  pass  through  the  Ej,  the  indices 
IPIVj  and  ij  are  generated  recursively  as  follows: 

for  j = 1 until  abs(NEL)  do 
begin 

i = HL(JL+j); 

if  i > 0 then  ( IPIVj  = IPIVj. i*.  ij  = i ) 

else  ( IPIVj  = abs(i);  ij  = IPIVj.,  ); 

end. 

For  FTRANL  the  indices  IPIVj  and  ij  are  used  as  soon  as  they  are 
generated.  However,  for  BTRANL  the  indices  IPIVj  must  be  stored  in  a work 
vector  during  the  above  forward  pass,  and  then  used  in  reverse  order. 

Given  from  point  (a)  above  that  NEL  is  typically  quite  small,  the  double 
handling  of  the  indices  HL(JL+j)  does  not  add  up  to  a significant  overhead, 
and  in  a pure-Fortran  implementation  the  overall  storage  management  is 
simplified  by  retaining  just  the  two  parallel  arrays  HL,DL.  (In  a machine- 
coded  implementation  the  efficiency  of  FTRANL  and  BTRANL  would  be  readily 
improved  if  the  IPIVj  for  all  updates  were  stored  permanently,  and  the 
increase  in  storage  would  not  be  substantial.) 


i 

i 
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6.2.2  Storage  for  U 


REALMS 

DU(MAXU) 

Packed  nonzeros  for  each  column  of 

U. 

INTEGER*2 

HU(MAXU) 

Corresponding  row  indices. 

INTEGER 

NU(M2) 

Number  of  ncnzeros  in  each  column 

of  U. 

INTEGER*2 

HPIVL(M2) 

Pivot  rows. 

Related  INTEGERS: 

MAXU 

Maximum  number  of  nonzeros  currently  allowed  for 
U file.  (Reset  during  each  call  to  INVERT.) 

in  the 

M2 

Previously  defined. 

NETAU 

Current  number  of  columns  in  U,  including  those 
have  been  deleted. 

that 

KDU,  KHU,  KU 

See  also. 

The  N-th  column  of  U is  represented  by  NU(K),  HP!VU(K)  and  a set  of 

contiguous  entries  HL'(JU+j),  DL;{JU+j)  •f'or  j = 1,2 NU(K).  Recall  from  the 

definition  of  U and  F that  each  column  of  U represents  a conventional  eta 
transformation  with  pivot  element  1.0  (which  is  therefore  not  stored). 

When  U is  updated,  new  columns  are  added  to  the  U file  in  the  same 
format.  If  a spike  was  deleted  during  the  basis  change,  the  corresponding 
column  0^  II  is  flagged  by  changing  the  sign  of  the  entry  in  HPIVU.  The 
storage  for  this  column  is  not  recovered.  If  a triangle  column  was  deleted 
from  the  basis,  any  packed  nonzeros  in  the  relevant  pivot  row  are  physically 
reset  to  zero. 
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6.2.3  Storaae  for  F 


The  elements  of  the  upper  triangular  matrix  F, 
stored  column-wise. 


REAU8  F(NF) 


Related  INTEGERS: 


MSPK  Maximum  number  of  columns  currently  allowed  for  in  F. 

(Reset  dulling  each  call  to  INVERT.) 

NF  Dimension  of  array  F,  namely  MSPK*(MSPK+l)/2. 

NSPIKE  The  number  of  spikes  found  by  P4.  Later,  the  current 
number  of  columns  in  F (set  in  FACTOR  and  MODLU). 

KFF  See  also. 


If  NSPIKE  = 3 and 


F = 


2 5 1 

8 0 
9 

.j 


then  array  F contains  3*(3+l)/2  = 6 elements,  namely  2,  5,  8,  1,  0,  9 
in  that  order. 

If  a spike  is  deleted  during  a basis  change,  the  corresponding  row  and 
column  of  F are  deleted  and  the  storage  is  recovered. 
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6 . 3 Cholesky  Factor  of  the  Reduced  Hessian 
REAL*8  R(NR) 


The  elements  of  the  upper  triangular  matrix  R, 
stored  column-wise. 


Related  INTEGERS: 


MAXS 

MAXR 


NR 


The  maximum  number  of  superbasic  variables  allowed. 

The  dimension  of  the  reduced  Hessian  factor,  R. 

1<MAXR<MAXS.  If  the  number  of  superbasic  variables 
exceeds  MAXR,  the  information  in  R is  no  longer  updated. 

An  exception  is  in  DELCOL;  if  the  JQ-th  column  is  being 
deleted  and  JQ<MAXR,  then  R is  updated  normally,  although 
R will  not  be  used  unless  the  number  of  superbasics  decreases 
to  MAXR.  (It  vwuld  probably  be  preferable  to  enforce  R = I 
whenever  the  latter  occurs,  and  skip  any  earlier  updating  in 
DELCOL.) 

The  dimension  of  array  R,  namely  MAXR*(MAXR+l)/2. 


Since  R is  in  general  a dense  triangular  matrix,  it  is  stored  in 
exactly  the  same  manner  as  the  matrix  F on  the  previous  page. 


7.  INVERT  STATISTICS 


The  following  items  are  printed  in  the  Iteration  Log  whenever  the  basis 
matrix  B is  re-factorized.  The  factors  L,  U and  F have  been  defined  in 
section  6.2.  The  constraint  matrix  array  is  denoted  by  A. 


Label 

FACTORIZE 

DEMAND 

ITERATION 

INFEAS 

OBJECTV 


SLACKS 

LINEAR 

NONLINEAR 
ELEMS 
DENSITY 
P4  BUMPS 

SPIKES 

CORE  REQD 


Description 

Counts  the  number  of  factorizations  since  the  start  of 
the  run. 

Gives  the  reason  for  the  present  factorization,  as 
described  in  section  8.  (Variable  INVRQ). 

The  current  iteration  number. 

The  number  of  infeasibilities  at  the  start  of  the 
previous  iteration. 

If  INFEAS  > 0,  this  is  the  sum  of  infeasibilities  at  the 
start  of  the  previous  iteration. 

If  INFEAS  = 0,  this  is  the  value  of  the  correct  objective 
function  after  the  previous  iteration. 

Ine  number  of  slack  variables  in  the  basis. 

The  number  of  basic  structural  variables  that  are  not 
involved  in  the  nonlinear  objective. 

The  number  of  nonlinear  variables  currently  in  the  basis. 

The  number  of  nonzero  matrix  coefficients  in  B. 

Percentage  nonzero  density,  100*ELEMS/(M*M) . 

The  number  of  bumps  or  blocks  found  by  P4  while  permuting 
B to  block  lower  triangular  form. 

The  number  of  spikes  tentatively  assigned  after  applying 
P3  to  each  bump  (prior  to  LU  factorization). 

The  number  of  words  of  core  required  for  P4  and  P3  to 
determine  the  bump/spike  structure  of  B.  This  will 
always  be  substantially  less  than  the  total  core  required 
for  the  run  (since  the  subsequent  LU  factorization 
requires  significantly  more  core). 
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L LIMIT 


1 


I; 

ij 

I 

If 
f- 


U LIMIT 
LU  BUMPS 

SPIKES 

AIJ  ELEMS 


L ELEMS 

U ELEMS 

F ELEMS 


The  number  of  nonzeros  allowed  for  in  the  LU  factoriza- 
tion and  subsequent  updates  to  L. 

Similarly  for  U. 

The  number  of  bumps  after  LU  factorization.  This  will 
invariably  be  the  same  as  P4  BUMPS. 

The  number  of  spikes  after  LU  factorization.  Column 
interchanges  may  increase  the  number  tentatively  assigned 
by  P3  above. 

The  number  of  nonzeros  in  L that  are  represented  by 
pointers  into  A. 

These  belong  to  non-spike  (i.e.  triangle)  columns  of  the 
basis,  and  will  often  comprise  the  bulk  of  L,  unless 
there  are  very  many  spikes.  This  mode  of  storage  can  be 
suppressed  by  a card  reading 
IMBED  NO 

in  the  SPECS  file.  The  L ELEMS  count  will  then  corres- 
pondingly increase. 

The  number  of  nonzeros  stored  in  the  L file.  This 
includes  one  for  each  column  that  is  imbedded  in  A as 
above. 

The  number  of  nonzeros  stored  as  packed  column  vectors 
in  the  U file. 

The  number  of  nonzeros  in  the  (dense)  triangular  matrix 
F,  followed  by  the  percentage  density  of  F.  These 
figures  are  for  interest  only;  they  indicate  the 
efficiency  or  otherwise  of  treating  F as  a dense  matrix. 


The  following  statistics  refer  to  the  LU  factorization  of  B.  They  are 
printed  only  if  one  of  the  items  TRISWAPS,  SPKSWP  or  REJECTED  is  nonzero. 

TRISWAPS  The  number  of  times  a triangle  column  was  interchanged 

with  a spike  in  the  same  bump  to  preserve  numerical 
stability.  The  LU  SPIKES  count  usually  exceeds  P4  SPIKES 
by  this  amount. 
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SPKSWP 


REJECTED 


MIN  PIV  RATIO 


TOLS 





The  number  of  times  two  spikes  in  the  same  bump  were 
interchanged  to  preserve  numerical  stability. 


The  number  of  variables  rejected  from  the  basis  and 
replaced  by  suitable  slacks.  (This  seldom  occurs,  unless 
the  basis  supplied  at  the  start  of  a run  proves  to  be 
singular,  or  else  the  current  B is  very  ill-conditioned.) 
Any  such  rejected  variables,  up  to  a maximum  of  20,  are 
stored  in  a list  and  given  priority  for  reintroduction  to 
the  basis  during  subsequent  calls  to  PRICE. 

The  smallest  preassigned  pivot  ratio  encountered  and 
accepted.  (This  will  always  exceed  LU  ROW  TOL.) 

For  each  bump  column,  the  pivot  ratio  is  measured  as  the 
size  of  the  preassigned  pivot  element  in  the  correspon- 
ding column  of  L,  relative  to  the  largest  element  in  the 
same  column  of  L and  inside  the  relevant  bump  (including 
the  pivot  element  itself). 


The  values  just  used  for  LU  ROW  TOL  and  LU  COL  TOL 
respectively.  (If  numerical  error  is  apparent  following 
LU  factorization,  the  factorization  will  be  repeated 
with  larger  values  for  these  tolerances.) 
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8.  INVERT  REQUESTS 


The  variable  INVRQ  is  set  in  DRIVER  and  certain  other  subroutines  to 
request  a call  to  INVERT.  If  two  such  requests  occur  without  an  iteration 
being  performed,  an  error  exit  is  taken  from  DRIVER.  If  INVRQ  > 20  an 


iteration 

was  beginning  but  could  not  be  completed. 

INVRQ 

Set  In 

Meaning 

0 

DRIVER 

Start  of  run. 

1 

DRIVER 

Invert  frequency  interrupt. 

2 

•3 

DRIVER 

LU  file  has  grown  significantly. 

o 

4 

MODLU 

No  room  to  update  L. 

5 

MODLU 

No  room  to  update  U. 

6 

MODLU 

No  room  to  update  F. 

7 

8 

MODLU 

New  pivot  element  in  F is  too  small 

9 

10 

DRIVER 

Row  error  after  row-check  interrupt 

11 

12 

13 

14 

15 
21 
22 


DRIVER 

RGITN 


No  workspace  in  LU  file  for  use  by  CHUZR. 
Linesearch  failure. 


9.  ERROR  CONDITIONS 


The  variable  lERR  is  set  by  many  subroutines  to  indicate  an  error 
condition.  Some  values  cause  the  run  to  terminate.  If  lERR  < 10  it  is  safe 
to  output  the  final  basis  and/or  print  the  final  solution. 

I ERR  Set  in  Meaning  Suggested  Remedy 

-4 


-1 

L^ITN 

With  composite  objective,  the 
incoming  variable  will  not 
improve  the  true  objective. 

0 

DRIVER 

Proceeding,  or  solution  is 
optimal . 

1 

DRIVER 

Problem  is  infeasible. 

Check  constraint  data. 
Raise  FEASIBILITY  TOL. 

2 

LPITN, 

RGITN 

Problem  is  unbounded. 

Add  realistic  bounds  to 
variables. 

3 

DRIVER 

Too  many  iterations. 

Continue  run. 

4 

DRIVER 

INVERT  called  twice  in  a row. 

System  error. 

5 

RGITN 

Too  many  superbasics. 

Increase  SUPERBASICS. 

6 

SETX 

Row  error. 

Scale  data  better. 

7 

RGITN 

Unable  to  improve  objective 
value  after  repeated  attempts. 

Is  the  nonlinear  objective 
and  gradient  well-defined? 

8 

DRIVER, 

FUNGRD 

User's  CALCFG  •requests  that  the 
run  be  terminated. 

9 

RGITN 

CHUZQ  failed  too  often  to  find  a Check  scaling  of  the  data, 
superbasic  to  replace  a variable 
that  wants  to  leave  the  basis. 

10 

INVERT 

Not  enough  core  for  P4.  j 

11 

12 

FACTOR 

FACTOR 

Not  enough  core  for  L. 

Not  enough  core  for  U. 

Increase  core  by  a 
substantial  amount. 

13 

FACTOR 

Not  enough  core  for  F.  i 

14 

INVERT 

Wrong  number  of  basic  variables. 

System  error. 

15 

16 

17 

18 

DRIVER 

Not  enough  core  tor  INSERT. 

Increase  core. 

19 

DRIVER 

Not  enough  core  for  LOADN. 

Increase  core. 
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I ERR 

Set  in 

Meaning 

Suggested  Remedy 

20 

LOADB 

Incompatible  basis  map  supplied. 

Check  OLD  BASIS  FILE. 

21 

DRIVER 

Not  enough  core  to  output 
solution. 

Increase  core. 

22 

23 

24 

25 

In  general,  error  conditions  and  warnings  are  brought  to  the  user's 
attention  by  a printed  message  starting  with  XXX  . 

Caution: 

The  message  EXIT  --  OPTIMAL  SOLUTION  FOUND  must  always  be  interpreted 
in  terms  of  the  size  of  the  reduced-gradient  vector,  particularly  if  messages 
of  the  form 

XXX  NO  FUNCTION  DECREASE. 

occur  at  the  end  of  the  run.  The  final  computed  solution  may  or  may  not  be 
acceptably  close  to  optimal,  depending  on  the  "friendliness"  of  the  nonlinear 
objective  and  the  condition  of  the  current  basis.  The  latter  message  occurs 
in  RGITN  if  the  linesearch  procedure  was  unable  to  obtain  an  improved  objec- 
tive value,  at  a time  when  various  tolerances  suggest  that  it  should  have. 

The  following  information  is  also  printed: 

NORM  RG  = the  largest  reduced-gradient  component  amongst  superbasics. 

NORM  PI  = a measure  of  the  size  of  the  current  tt  vector. 

LAST  SB  = the  column  no.  of  the  last  superbasic  variable.  This  was 
previously  nonbasic  with  reduced-gradient  shown. 

The  recovery  action  taken  by  MINOS  is  to  add  the  most  favorable  nonbasic 
variable  to  the  superbasic  set  and  continue  ootimization  in  a larger  subspace. 
If  PRICE  finds  there  are  no  favorable  nonbasics,  the  solution  is  accepted  as 
the  best  attainable  and  is  declared  OPTIMAL. 

Broadly  speaking,  if  (NORM  RG)/(N0RM  PI)  = 10"*^,  then  the  objective 
function  would  probably  change  in  the  d-th  significant  figure  if  optimization 
could  be  continued.  The  user  must  decide  for  himsel f whether  d is 
sufficiently  large. 


10.  SUBROUTINE  SPECIFICATIONS 


Each  of  the  subroutines  in  MINOS  is  documented  under  the  following 
headi ngs : 

i 

• Subroutines  called  and  called  by 

• Purpose 

• Parameters 

■ • Method 

The  Parameters  description  does  not  include  any  variables  with  the 
same  name  and  function  as  those  defined  in  section  4. 
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SOBFOOTIHS  ADCCCl  ( 2,MNN,IIJ:,NP  1,  NF,  NP , KN , KS, K S,  f! AXt, f5 A Xtl , 

1 HSPIKE,HA,PE,  PPIVL,HP1  V'),HPI  VF  ,HE,KL,  Ni;, 

2 A,F,X,GPD,T,  Y 1,XN,C,G1,GNEH,HL..Ol,F!(l,  C'l,F,iINIT,  JQ  ) 


Subroutines  calletl:  Called  ly: 

PTRABO  CALCG  PTPANL  FUNGPD  EGITP 

GETGED  SETPI  OKPACK 


Purpose: 

To  construct  a new  coluinn  for  F when  th«  no,  of  superfcasics  is  hoir.q 
increased  by  1,  froir,  NS-I  to  KS. 


Faraneters: 

ISTEGEF  JQ  (Input)  Cclumn  no.  of  variable  being  add'-^d  to 

the  set  of  supertasics. 

LOGICAL  OMIT  (Input)  Specifies  if  new  column  of  R is  to  be 

a unit  v«c»-or  or  net. 


Method : 

If  OBIT  = TPOE,  the  unit  vector  c (NS)  will  be  used.  This  is  tie  case 
if  BUITIPIS  PPICEing  is  inveked. 

If  ONII  “ FALSE,  an  atteapt  will  be  aade  to  estieate  curvature  in  the 
new  sutspace  by  f inite-dif fersrclng  the  gradient  v^sctor  alona  the  new 
colufn  of  Z (the  null  space  na'-rix).  See  flurtaoh  and  Saunders  ( 1978  ), 
section  3.1.4. 

If  the  off-diagonal  elewents  cf  the  n^w  cclurn  show  signs  of  being 
extreeely  big  or  stall,  fur«:hai  work  is  avoided  and  th'*  unit  v'*ctor 
e (NS)  is  used  instead.  If  the  rev  liaqcnal  eieirent  proves  to  be 
significantly  larger  cr  sealler  than  the  ncr®  of  F,  again  a urit 
vector  is  used  and  the  werk  required  ♦:o  estinate  the  new  coluan  is 
wasted.  The  aim  here  is  to  avoid  cons*"  ructing  a n^w  F <■  hat  is 
significantly  wore  ill-conditioned  than  the  previous  R. 
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SrJBBOOTlME  ALIGN!  N , NS  , HE  , BL,  3L’ , X,  T , IHFr  A3  , MALIGN,  1 SPE^  L I 


Subroutines  called:  Called  hy: 

DOT1  RGITN 

Purpose : 

To  scale  the  search  iirectioE  for  superbasics  so  that  if  sev»^Lal  are 
about  to  reach  a bound  then  all  of  them  do  so  siiultaneous  ly. 


Paraaeters: 


INTEGEP*2 

INTEGFP 

BEAl*3 

BFAI*8 

REAI*B 

SEAl*8 


8B(1»S) 

malign 

THETAS  (NS) 
TSPBCL 
X (MS) 

Y (MS) 


(Input)  Cclujir  nos.  of  superbasics  (orlv). 
(Output)  Nc.  of  variables  finally  aliened, 
Workspace  to  store  stens  to  reach  bounds. 
(Output)  Firal  step  si??,  if  NALION>0. 
(Input)  Values  of  superlasics  (cily). 
(Input)  Search  direction  for  superbasics. 


Kethcd: 

Classify  each  supertasic  variable  X (j)  as  "norinal"  or  "sT^ecial" 
as  follcus.  If  the  seitch  direction  Y (1)  is  larger  than  scti'a 
tolerance  ZTCLY,  and  if  the  distance  *0  the  relevant  bound 
(depending  on  the  sign  cf  Y(i))  is  less  than  tcme  tolerance  ZTCI.'', 
then  X(j)  is  special,  otherwise  noraal. 

2.  Ccnpute  TNOBHt  (or  I3FECI),  the  step  that  lakes  thn  first 

noraal  (or  special)  variable  hit  its  bound.  F^ore  cteps  for  all 
special  variables  ir.  array  ItFlAS. 

3.  If  the  first  supertasic  to  reach  a bound  is  special  (i.e.,  if 
TSPFCL  <=  TMCPfL)  , scale  dctin  the  coioonents  Y (j)  for  special 
variables  so  that  all  special  variables  will  reach  thair  bounds  if 
a step  TSPECL  is  taken  along  the  icdified  search  direction.  (A 
basic  variable  tay  reach  a fourd  earlier,  in  which  case  the 
alignient  will  have  been  tc  no  avail.) 


SOBBOUTIHE  BTRANL(  *1  , M 2 , N E , HP  1 , H AXL  , HA  , h h,  A , H PI  VL , HL  , HL  , PI  , HW  , Y > 


Subroutines  called:  Called  L7: 

Ncne  CHU7.0  FACTOP  PRICE 


Purpose: 

To  solve  the  systea  1 (ttan'ipose)  *x  = y,  Mbf-rT-  both  » •nd  y are 
stored  in  the  paraaeter  Y,  This  is  t.he  2r.d  cart  of  r rcn vt  ntior.al 
"ETHAN"  operation. 


Faraaeters: 

r TTEGE8e2  HW (H)  A work  vector  required  to  re-qenerate  the 

sequence  of  row  interchanqes  that  occurred 
dorinq  each  update  tr  the  LU  factors.  (Hot 
required  during  INV.EPT)  . 


Kethod : 

The  reverse  of  FTBANl.  Updates  are  treated  first,  followed  by  *he 
t ransforaations  constructed  by  FACTOB. 

Treataent  of  the  update  transfer aar ions  is  rather  complex  because  just 
o;.e  integer  is  stored  per  nonzero,  rather  than  twe,  and  a forward  pass 
through  the  integers  ir.  each  ttansf  or  nation  is  reu’iired  *^0  re  lenerate 
the  second  set. 

Puring  INVERT,  otly  a partial  sweep  through  the  noraal  colupr. 
transformations  is  needed.  The  variable  LFGN  allows  the  backward 
sweep  tc  stop  at  the  beginning  of  che  current  bump.  It  is  reset  ♦■o 
point  to  the  beginning  of  1 before  exit  from  FACTOR. 


SUBPCUTIIiB  BTSHNOI  H,  K 2,  NF,  H AX  (J,  HFI V tl,  JT  TV  !•  , NU,  HO,  DU,F  , Y ) 


Sabroutiass  called:  Called  by: 

None  ADCCll  LPITH  RGITI  SFTX 


Poryose: 

To  solve  the  system  l)*i  = y,  where  both  x and  y ace  stocf  1 in  the 
j paraietec  Y.  This  is  the  2nd  part  of  a conventional  "FTBAN" 

i operation. 


Hat  hod: 

Normal  column-wise  back-substitution,  with  each  column  of  U and  F 
being  treated  as  a column  of  the  full  matrix  ”0".  Columns  are 
processed  backwards.  If  a ccmpu^-od  cowponent  of  x is  ''ssratially 
2ero,  the  corresponding  coluins  cf  U and  h can  be  skipped. 
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SUBFCUTIM2  BUMPS  ( li  , N 1 , N Z , IF  , IBN  , CUTF,  NUMB,  3,L0HL,  PBEV  , NBLK  ) 


Sobroutioes  call6<1:  Called  by: 

None  Fu 


Put  pose : 

Given  tha  roe  nuebers  ol  the  ncnjctos  in  each  column  of  a square 
sparse  matrix  A,  this  routine  finds  a permutation  U such  ♦hat  ?*Ai' 
is  tlock  upper  triangular. 


Parameters: 


INTEGEF 

N 

(Input)  Dimension  of  the  matrix. 

ISTEGBP 

N1 

(Input)  N*1. 

INTEGE? 

NZ 

(Input!  Bo.  of  ncnzeros  in  the  matrix. 

INTEGES*2 

IP(N1) 

(Input)  Feints  +o  the  first  nonzero  in  each 
column. 

IVT5:r.EP*2 

IRN  (N7) 

(Input)  List  of  row  numbers  for  nonzeros. 

I'»TEGE?*2 

OOTF (N) 

Work  vector. 

INTEGBP*2 

NOBF  (N) 

(Output)  NUMD(k)  will  contain  the  permiitfed 
pcsition  of  column  K. 

r.«TEGEP*2 

E(R> 

(Output)  B(I)  will  contain  the  row  no.  in  ♦he 
permuted  matrix  of  the  beoinninq  of  ♦he  I~th 
block  . 

: lTtGEB*2 

LOHL(N) 

Work  vector. 

; .•"2GEF*2 

FBBVL  (N) 

werk  vector. 

; 1 “EGEF 

NBLK 

(Output)  Bo.  of  blocks  found  (includinp  ores  of 
dimension  1). 

f ■ thod: 

£ne  C.  Tarjan,  SIAM  J. 

Ccxputirg  (1972),  1,  pp.  146-160. 

£ ibrcutine  BOHPS  is  derived  directly  from  Harwell  subroutine  ''CllA  bv 
C . , Peid  (1975).  (MC13A  has  sinc^  heor.  replaced  in  th«  Harwell 
f Jbroutine  Library  ty  MCHC;  see  l.S.  Duff  and  0.  K.  Feid,  "An 
i onljnentation  of  Tarlan's  algorithm  for  the  Mock  trianaul arization 
( a matrix.  Report  No.  CS5  29,  Harwell,  April  1976). 

It  is  assumed  that  the  matrix  A has  already  been  permuted  to 
hove  ncnzeros  or.  its  diagonal.  Sec  the  transversal  finder,  sihroutine 
TPNSVL. 
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SUEFCtJTIHE  CALCFG(  WODE  , > ,X,P  ,G , NST  ATE , N E ROP  ) 

Sobcoutines  called:  Called  by: 

Hone  FiiNGPn 


Parpose: 

To  evaluate  a user's  particular  nonlinear  objective  function  f (x) 
its  qradiect  vector  g(x)  at  a given  point  x. 


Paraireters; 

The  required  fora  for  this  sutroutine  is  docueented  in  tiie 
Oser's  Guide,  section  III. 3. 


and 


MINOS 
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SUBBCUTIBE  CMCG  ( H , H F ,NP  1 ,IIS  , 1*3  , , HB  , H?  , X , ?!  ,>3BD,G,  BOSOI  K ) 


3ubrcutines  called:  Called  by: 

None  XDDCCL  OBIVFF  FGITN 

Purpose: 

To  calculate  the  reduced  gradient  vector  G for  tli-a  current  t of 

superbasic  varlableF. 

Parameters: 

5PAL*8  GRD(nS)  (Input)  Ccntains  the  current  gradient  vector  for 

superbasics  in  positions  ...»  M<-N5. 

FPAL*8  GGNORH  (Output)  Feturns  the  norT  of  the  reduced 

gradient  vector.  This  Is  taken  to  he  *he 
largest  elertert  cf  G in  absolute  size. 

het  hod : 

The  J-th  ccsponent  cf  the  reduced  gradient  is 
G(J)  = GPD(?!^J)  - PT*A(J) 

V lere  A (J)  is  the  cotrespccding  cclunn  of  A. 


l; 


SHBFOOTINS  CHKCIF(  K , BL , BU , X , C , H E EA S ) 


Sobicutines  called:  Called  by: 

Hone  CHKG!^r 


Per^cse ; 

CHKGPO  is  about  to  use  direction  D for  forward  dif ferer.rirq  'he 
nonlinear  ohjectiwe  functicn,  starting  frow  the  first  feasi  i le  pcir.t 
X.  This  routine  checks  that  D is  a feasible  direction,  and  ..  titles  D 
if  necessary. 


Paramaters: 


IHTEG?F  N (Input)  Will  be  KN,  the  no.  of  nonlin'=ar  vars. 

BL(N)  (Input)  Lower  fccunds  for  the  ncnlin^a'  vars. 

BO(N)  (Input)  Upper  bounds  for  the  nonlinear  vars. 

PEAL*B  X(H)  (Input)  Will  be  TN,  the  nonlinear  variibles. 

HERL^rt  E (H)  (In, out)  The  diroctior,  to  he  checked. 

TNTEG3P  NPEAS  (Output)  The  no.  of  feasible  c^iirponeut.';  in  the 

(possibly  BCdified)  vrctcr  D. 


SOEFOOTIME  CHKGPD*  N ^ BL , BiJ  , . G , T A 1 1.  , ft  . H , C , FPSHCH  ) 


Subroutines  called:  by: 

CALCFG  CHKDIF  DOT  rFIVFF 


Purpose : 

To  verify  to  soae  extent  that  the  user  has  rerrectly  programnitd  the 
gradient  of  the  ofclective  function  in  suhroutine  CALCFG. 


Parameters: 

INTEGER  N (Input)  Will  he  NN,  the  no.  of  nonlinear  vars. 

HEAL  BL(N)  (Input)  Lever  bounds  for  the  nonlinsdr  vars, 

FEAL  bo  (N)  (’’Tjput)  Upper  fcctinds  for  the  nonlinear  vars. 

REALMS  X (N)  (Input)  Will  be  fN,  the  r.onlinedr  variables. 

REAI*8  G (S)  (Input)  Will  be  C,  th“  gradient  at  th--  first 

feasible  point. 

INTEGPR  IP&IL  (Output)  Will  be  zero  if  the  gradient::  lock  CK, 

2 if  they  don't,  or  B if  the  us'^'-r  sets  ♦o 

a negative  number  in  subroutine  CAl'FG,  to 
teroinate  the  run. 

FEAL*8  A (N)  fcork  vector  for  1st  differenca  directicn, 

HEAL*8  B (N)  Work  vector  for  2nd  difference  direction. 

FEAL*8  C(ll)  Work  vector  to  held  gradients  at  XtH^A, 

Contents  not  used. 

I £AL*8  EPSKCH  (Input)  Machine  precision. 


f jt  hed : 

This  routine  was  taken  directly  from  subroutine  C^’KGPD  l,y  Gill, 
Murray,  Picker  and  Barter.  S*-:e  rccuaant  no.  rU/05/0/F/1 1/  ^5,  DN^J, 
hational  Physical  Laboratory. 

S'iuht  changes  were  made,  e.g.  tc  make  sure  A and  B were  fed'sibi*' 
c .r''cticns. 
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SU3F0UTIKE  CG  ( K S , G , G N F W , r , T I:  ET A ) 

Subroutines  called:  Called  by; 

DCT  DOT1  FGITN 


Purpose; 


To  compute  a search  directicn  for  supei  tas  ics,  usinq  c r.c-  cf  severcl 
conjugato  gradient  mctf-ods  for  unccnst  rai  ne  ! mlnitnizat  ion. 


Parameters: 


REAL*8 

PEAL*8 

?EAL*8 


G (US) 
GNEW  (NS1 
C (NS) 


REAL*8 


THETA 


(Input)  Seduced  gradient  at  previous  s+ap. 
(Input)  Eeduced  gradient  at  current  s''-^d. 

Cr.  input,  the  previous  search  direction  for 
Euperlasics. 

Cn  output,  tli.3  nt  w search  direction  for 
Euperkasics. 

(Input)  The  step  moved  during  the  prt  vious  itn 
along  the  (old)  search  direction  D. 


Method; 

1.  Use  MSGCG  to  print  a message  if  this  is  the  first  entry  t^.  ^3, 

2.  Use  ITNCG  to  determine  if  NS  consecu  ive  iterations  hive 
performed  since  the  last  restart.  If  so,  request  a restar*. 

Use  MOnCG  to  determine  which  version  of  CG  has  been  requested. 

h.  Host  methods  (Fletcher-Raeves,  Polak- Pibiere,  Perry)  compute  th ^ 
new  D according  to 

r (J)  = -GNEW(J)  ♦ CG2tTA*D(J) 

for  some  scalar  CGEETA.  Heioryless  CPF,  CC"DFP,  etc.  add  a terra 
GAHnA*  (GNEW  (J)  - G (J)  ) 
for  some  scalar  GAMKA. 
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il 

j( 


Note;  CG  is  called  from  RGITN  at  tha  end  of  an  iteration,  bu*  only  if 
that  iteration  was  unccnst raine d . Tha  required  search  directicr.  D is 
then  still  available  for  us*^  during  the  next  call  to  FGITN.  '■o'-e  that 
in  PGITN,  n resides  in  Y(Mel),  ...»  Y(K*NS).  This  is  overwri-ttn  if 
the  ensuing  step  is  constrained. 
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SOBFOOTIN’!  CHU7.Q  ( P,  N , SE,  NFI,  MS, 

1 HA,HE,HE,A,El,BO, X, Y,JO,riVCT,VAF HE!  ) 


':«brcutines  called:  Called  bv: 

FTPANfJ  BTBAKl  FfllTK 


Purpose . 

When  the  JP-th  variable  leaves  the  basis,  thi§  routine  selects  one  of 
the  supertasic  variables  to  rejlac?  it.  The  prirnary  aim  is  to  k»ep 
the  basis  well-conditicred.  where  possihlr  the  choice  is  tiase’ 
towards  a superbasic  that  is  away  ftoa  its  bounds,  in  the  hot's  that  it 
will  be  less  likely  to  leave  the  basis  durinq  cjbsequent  iterations. 


Faraveters; 

Y (MS)  The  first  b components  of  Y contain  * Ue  JP-th 

tCH  cf  B inverse.  On  output,  Y(I'+1),  .... 

y JM+NS)  will  certain  a vector  v reauired 

later  by  BIPBOb, 

JQ  (Output)  Will  pcint  to  the  chosen  superbasic. 

b < JC  <=  MVMS, 

PIVOT  (Output)  The  value  of  Y**J.  (K)  for  th.j  chosen 

cclunn,  K = FE  (JQ) . 

VARMET  (Input)  Used  tc  skip  coirputing  v if  ir  ic  no* 

goitej  tc  be  used  later. 

■ f. 


; . * hed  ; 


As  described  in  Murtaqh  and  Saunders  (1978),  section  1.3.2. 


REMV8 

INTEGER 

REAL*a 

:OICAL 


SnnROUTINE  CHDZ?  ( »',(i,HP,PL,EI!,X,Y,  ATCN'f^,nJO,TOLy,TOLPTV,'"f'Z,  a,J 
1 RTFETA  ,TAL' ,I  F''Ay  1 


Subroutines  called:  Cilltd  ^y: 

None  IPITN 


To  choose  an  indax  JP  and  a step  siz*=  TFEIA  >-  0.0,  such  that  variahl* 
HE(JE)  reaches  one  cf  its  fc.ciirds  when  X chanqes  to  X - THETA*y.  If  K 
is  feasible,  H3  (JP)  will  be  the  first  varia)-le  tc  reach  a hour.  ;.  if  X 
is  infeasible,  JP  is  chesen  tc  ainiPize  the  su?  of  infeas i h i 1 i t ies 
without  reqard  to  the  numter  (which  may  ir.crtase). 


Paraaeters: 
INTEGER  R 


RF.AL*B  Y (R) 

lOJ.ICAl  ATBND 


RFAltB  DJQ 


5EALt«  TCLX 

TCLPIV 

PEAL*H  TPHTA 

INTEGEP  JP 


(Input)  This  will  really  ce  Jt1  when  CHL’Z?  is 
called  frcB  LPIIN.  (The  variable  e>T;t.:  ;iq  tt  <= 
hasis  will  be  in  [osition  K*’,  and  Yd  + D will 
te  1.0.)  It  will  be  K+ iS  when  call<-ci  from 
SGITM,  and  th“  last  N5  values  of  Y will  (•  th® 
search  direction  fci  su perb-isics. 

(Tnpit)  "^he  search  direction. 

(Input)  IBIJE  for  LPITh  if  variable  er.tarir.u 
hasis  is  cttrertly  at  i*s  upper  bturi".  Ir, 
’ffeef,  Aif-sr  causes  -Y  to  1“^  used  in  ;iac'>  f 
Y in  the  usual  ratio  test.  In  qenetai, 
shculd  therefore  be  TP'JE  if  is  bainq  c'-^nq-') 
tc  X ♦ ThETA*Y. 

(Incut)  This  is  used  only  if  the  curie::*  X 's 
infeasible.  It  is  the  slop-  of  the  pi» c -wise 
linear  qraph  cf  the  sum  of  in  f€asit)il  iti- s 
versus  ThETA.  Rrr?  generally,  it  sh~il(J  h.' 
q'p,  the  usual  inner  pictuct  of  the  '*  ctive 
qradient  witii  the  search  direction,  p~  i . 

Foe  LFIIN,  this  is  the  reduced  gradien*  of  th^ 
inccffiinq  variable. 

For  FGIIN,  it  is 

q'p  = g*  (Zq) 

= fo'2)*q 

= (g'Z)t(-2'g) 

= - II  7' o It  (squared) 

= ~ ( I G II  (squared) 

whert>  G ccntair.s  the  redur*  3 qradient,  Z'g; 
q is  '■he  Dart  cf  p ‘•'or  su  perl  »s  i.-;  , 

namely  th5  negativ  rfducad  ■-l  idier.*; 
and  7 is  the  null-space  matrix. 

(Input)  Feasibility  tolerance. 

(Input)  A component  of  Y will  be  tr-i*-'d  r.s 
Zfctc  if  smallar  than  this  ♦<  lerance. 

(Cutput)  Th^  chesen  step  size. 

(Output)  Th-  ch'sen  iriax. 


IKTBG2P*2 


IKTI'GBP*2 

PFM*8 

BEAL*8 

IKTEGEP 


BPT(IFMAX)  (Wcrkspacc)  llcldj  foin’-ars  to  maintairi  a linked 
lif.t  in  th<"  attays  “RT,  FTHFTA  and 

TAD.  The  list  i’"-  oi  lertd  according  *0 
increasing  valuas  of  ‘At),  starting  frjB  ♦ha 
nagative  val'i*:,  ~a^o(l.)0). 

HPT  (IFMAX)  (Horurpvica)  I’oid?  nivct  cows  corresponding  *-n 

the  stops  in 
RTHETA  (IFr AX) 

(Workspace)  Holds  list  cf  s»ep-sizes,  THETA. 
TAU(IFHAX)  (Workspace)  Hcldt  lis^  of  slopes.  Tn'(I)  is 

the  slope  after  a ster  RTI!i'"A  (I)  is  tdk<’r.. 
lEMAX  (Input)  The  laxiBiuir  length  allcwad  '’or  tho 

lirked  list.  storage  for  the  4 above  arcavs 
ceres  fren  unused  space  in  tlie  LU  fiin  and 
therefore  varies  [..■'tween  calls  to  GHtiZP.  It 
will  .jlways  be  at  l«?.st  ''  (otherwis-  INVF.S’^ 
will  ) 6 call c d) . 


Set  hod : 

The  following  two  methods,  as  implemented  by  John  Tomlin,  November 
197b: 

♦ If  feasible,  Paula  Harris's  two-pass  method  for  taking 
advantage  of  neat  ties. 


* If  infeasible,  bennis  Ratick's  one-pass  method  for  minimizing 
the  sum  of  in  feasibilities. 


Harris 


The  aim  here  is  to  take  advantage  cf  the  likelihood  that  stveral 
variables  may  reach  a bcurd  almost  simultaneously  when  the  final  step 
TFFTA  is  taken.  (This  is  especially  probable  or.  degenerate  stens, 
i.e.,  when  TfiETA  = O.C)  The  method  does  net  expect  to  find  exr'ct  ti^s 
but  instead  looks  for  near  ties,  which  are  consid '•rably  more  frequent 
ai.d  mean  much  the  same  in  practice,  since  the  vertices  of  a -implex 
ar*-  in  seme  sense  "fuzzy",  at  least  tc  the  extent,  ct  the  uncer tain^ y 
in  the  original  data. 

1.  (First  pass)  Compute  a quantity  ESI  such  that  at  least  one  variable 

wculd  overshoot  its  bound  by  an  amount  FEPTHN  if  a step  rsi  were 
taken.  Since  FOPHC  regards  variables  as  being  feasibl-*.  only  ro 
within  a tolerance  TCLX,  we  won't  mind  if  certain  new  variables  go 
outside  their  bounds  by  a similar  amount.  PEFTBK  can  therefore, 

with  safety,  be  as  larga  as  TCLX.  Thinking  furthor  along  th-^se 
lines,  observe  that  the  actual  s-‘ep  THETA  that  is  finally  taken 

will  be  strictly  less  than  FSI.  Hence  it  is  OK  to  allow  PFr"'r>*;  to 

be  slightly  larger  than  TCIX.  we  take  the  value  .lore 

thinking  required  to  decide  what  the  maximum  value  should  really 
be.  (Actually,  the  real  reason  for  taking  PEFTdN  > TOLX  is  simply 
cn  principle,  to  avoid  getting  a zero  numerator  in  the  ratios 
defining  PSI.) 

2.  (Second  pass)  How  go  through  cemputirg  the  usual  ratios,  tliis  time 
witho'it  perturbation.  .skip  any  that  are  larger  ♦han  PSI.  of  the 
remainder,  find  the  one  for  which  the  pivot  elemer*-  Y (J>  is 
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larg<>st. 

Pacick 


A full  description  will  net  te  atteapted  here.  Briefly,  the  aethod 
would  to  conceptually  siapla  if  all  the  THPTA's  were  computed  and 
sorted  into  ascending  order.  Parick  showed  that  with  the  aid  of  a 
linked  list  the  sottirq  can  he  avoided  and  everything  can  be  doni^ 
during  a single  pass  thrciqh  X and  Y. 

For  further  details,  see  John  Tcnlin. 


I 
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SUBBOUTIRf:  CONOFPf  B , Y ,6  ,MB , N , IHETA  , ¥T  P,  RTP , KiDD, KSOB  ) 


Subroutines  called:  Called  by: 

H1MOD  BGITN 

Purpcse: 

To  ispleaent  the  "Coipleaentar y OPP"  or  "BPGS**  quasi-Newton  updat<>  to 
a given  Hessian  approxinaticn,  stored  in  factorized  fora  as  B'P,  thus: 

(new  F*R|  = (cld  8*B)  ♦ cl  yy»  ♦ c2  gg  * 

where 

g * new  reduced  gradient 
y » g - (old  g) 
p = search  direct icr 
cl  = 1/(THETA*y ’p) 
c2  = 1/g'p. 


Paraoeters: 

PEA1»8  « 

RPAl*8  I 

BEAL«a  I 

INTEGER  I 

INTEGER  I 


THETA 

(Input) 

The  step  just 

taken  along  direction  p 

TTP 

(Input) 

y p. 

GTP 

(Input) 

g*p. 

NADC 

(Output) 

Error  flag  for 

adding  the  rank-1 

r at li X, 

cl  yy*. 

KSUP 

(Output) 

Error  flag  fer 

subtracting  the  rank-1 

natrix. 

c2  qg*. 

Hethod : 

Call  RinOD  twice  with  appropriate  paraseters. 


: Subroutines  called:  Called  by:  jj 

fi 

f None  PRIVET  i' 

I . !i 

f 

Purpose: 

Tc  select  N vectors  froo  A in  a way  that  will  produce  a lower 
triangular,  nonsingular  basis  B.  Three  options  are  i»p lercsntcd , 
according  to  the  value  cf  IPAEM(I). 

IP9EM{1)  Reaning 

0 Set  up  the  all  slacK  basis. 

1 Scan  all  ccluans  cf  A. 

2 Scan  only  those  cclunns  of  A corresponding  to 
linear  variables. 


Paraieteis: 

INTEGEF*2  HS(N)  (Output)  Hill  define  chosen  basis. 

INTEGERV2  HPTVL  (M)  Horh  vector  to  record  state  of  rows. 


Method : 


1.  Initialize  state  vector  HS(-))  so  variables  are  nonbasic  at  their 
saallest  bound. 

2.  If  NS>0,  nahe  relevant  norlineat  variables  superbasic,  as  specified 
by  the  INITIAL  bcunds  set. 

1.  Look  at  the  values  given  to  nonlinear  variables  during  inpu^  of  the 
INITIAL  bounds  set.  Set  nonbasic  nonlinears  at  th^  bound  nearest 
to  the  values  in  XN,  ir  case  the  LO  or  UP  state  specified  during 
input  is  opposite  tc  the  default  in  1 above. 

4.  Proceed  to  construct  a triangular  basis.  The  array  HFIVL  will  have 
the  following  inter pretatlor : 


HPIVL  (i) 

Label 

Meaning 

K>0 

Pivoted 

Bov 

i has  been 

assigned  as  the 

pivot  tow 

for 

ccluin  K. 

-1 

narked 

Pow 

i occurred 

in  some  previous 

colusn 

but 

has  not  ye 

t been  assigned  a 

•s  a pivot 

row 

(i.e.  has 

not  been  labelled 

1 * voted) . 

0 

virgin 

Row 

1 has  net 

yet  occurred  in  a 

ny  chosen 

coluen. 
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Now  do  two  of  t.l)H  followlni),  flr«*  on  frmm  atroct 'ji  ^ j<i  only, 

th«n  on  all  *ti  net  ninlt;.  (lhi»r  (jotf  ti»  w cuctnr  alw  Into  tfi* 
hasiiw  flrat,  8iii>1«ct  of  coucoa  to  th«  tanlw  8r<iyinq  tciati'}"l«c.  f* 
is  quit*  I'ORKit  H rbAt  «oao  1 1 <j«  varxaiiitta  will  o^lll  ho  aotit.nufc 
Aftwi  CRASH.) 

(a)  For  any  given  ccluna,  find 

AIMAX  the  hi(jq*st  eleannt  in  tha  whole  coluan; 

AFIV(1)  the  hlciqe^t  elcBert  in  Racked  rows; 

APIV  (2)  thf  tiqqwst  alaaent  In  Virgin  rowt; 

NPIV  the  MO.  of  Pivoted  rows  in  this  coluan. 

The  rows  corrasrondinq  to  APIV(i)  are  the  "best"  Rarked  ai.d 
Virgin  cows  rasyertivelv. 

(t)  Select  a potential  pivot  cow.  If  possible,  take  the  best 

Virgin  row.  If  none,  and  if  NPTV  = 0,  choose  the  he^r  n^ckf-i 
row.  Otherwise,  skip  the  coluan  altogether. 

(c)  TOLPB  is  so»e  relative  pivot  tolerance,  e.g.  0.01.  To  prevent 
getting  a badly  conditioned  basis,  skip  the  col'itnn  if  the 
chosen  pivot  eleaant  (APIV(1)  or  ArTV(2)i  is  saaliet  than 
TCLPE*AI!1AX. 

(d)  Otherwise,  label  the  chosen  row  as  Pivoted,  set  the  state  of 
the  coluan  to  be  basic,  and  label  any  regaining  Virgir.  rows  in 
the  coluan  as  Marked. 

Finally,  fill  up  the  basis  with  the  Icgicals  on  any  rows  that  did 
not  get  labelled  as  Pivoted, 


SUBBOUTIME  DELCCl  ( M , NF , NS  , HS , HB , F , X ,JC , DEIST B ) 


Subroutines  called: 


Bone 


Called  ty; 

IBVEPT  PGITS  SQUEEE 


Purpose: 

To  perfora  housekeeping  when  the  JQ-th  superbasic  variable  is  deleted. 
This  aay  include  deleting  the  JQ-th  ccluan  of  P and  restoring  F to 
upper  triangular  fcra. 


Parameters: 
INTEGEF  JQ 


LOGICAL  DELETE 


Method : 


(Input)  The  position  of  the  suj erbasic  being 
deleted . 

(Input)  TEUE  if  quasi-Newton  or  coniugate- 
gradient  rethods  are  in  effect.  FALSE  if 
infeasible  or  ptcblem  is  linear. 


Perfera  a partial  forward  sweep  of  plane  rotations,  affecting  the  las* 
rows  and  columns  of  P.  Also,  move  ccluens  JC*1,  NS  of  R one 
place  to  the  left.  Adjust  arrays  H3  and  x similarly.  Decrease  NS  hv 

1. 
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SUBFOOTINE  DOT ( ) 

Subcoatines  called: 

None 


Called  by: 

CG  CHKGFD  3FABCH 


Purpose: 

To  coapute  tha  inner  product  £ 
V and  w. 


= V**w  for  the  tao  tJ-v«ctors 


Paraaeters : 

REA1*8  V(N{ 

FBAL*8  W(N) 

REAL*8  S 


(Input) 

(Input) 

(Output) 


S03BCUTINB  DOT1  ( N,V,S  ) 

Subroutines  called: 

None 


Called  by: 
ALIGN  CG 


FGITN 


Purpose: 

To  coupute  the  sub  of  squares  S = ((••¥  for  the  N-vector  V, 


Paraaeters: 

P!!AL*8  V(N) 

BEAL*8  S 


(Input) 

iCutput) 
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SOEBOOTINE  CPIVEF  ( Z , M MZ , IFSC LB , 

1 ^,H2,aS,ENh,B,BE,NE1,NF,NH,NM0,1AXS,NS, 

2 HSPIKE.r:E,HA,ft,El  ,PIJ,FT',XK, 

3 HS,HPIVl,HPIVO,HPIVF, 

4 BL,NU,U,C, X,tl, 

5 GPD,Y,T  1,G1,G,GNEH  ) 


Subroutines  called: 

Called  by: 

CALCG 

CHRGRD  CRASH 

DUHEN 

KINGS 

PCFSC 

FONGSD  GETGPC 

INSERT 

;•  1 

INVEST 

ITEROP  LOaCA 

ICACN 

1 ; 

LPI7N 

rCCLO  PRICE 

EL'NCH 

r 

RESETP 

RGITN  SAVES 

SETPI 

I’ 

SETX 

SOLN  STATE 

UNPACK 

r 

Purpose: 

i To  Invoke  all  necessary  algorithsic  routines  for  solving  a ptoblee, 

onca  it  has  been  input.  Ihis  includes  setting  up  a start  ina  basis, 
calling  INVERT  when  teauiied,  switching  between  simpiei  and 
; non-siaplex  iterations,  checking  that  tolerances  are  tight  enough  to 

declare  optiaality,  savinq  tasis  files,  and  cutputting  the  solution. 
Huch  cf  the  ccaplexity  cf  the  code  has  to  do  with  setting  the  Phase 
paraaeter  NPHS,  which  can  take  the  following  values: 


r- 

I 


; 


ij 

!; 

I 


I 


NPHS  Hea»’ing 


1 Ordinary  Phase  1 siaplex  aethod  (solution  iufeaalbl^) . 

Any  supertasics  will  be  ignored  by  PRICE. 

2 Ordinary  Phase  2 simplex  aethod  (solution  feasible) . 

Again,  any  superhasics  will  be  ignored  ty  PPICE. 

3 Ron-simplex  iteration,  (feasible  or  infeasibl*-)  , including 
a PSTCE  operation  to  select  additional  supertasics. 

4 Non-sinplex  iteration,  (teasikle  or  infeasible),  operating 
on  the  current  basics  and  supertasics,  without  a call  to 
PFICE. 

Osing  NXTPHS,  subroutine  POITN  suggests  a value  for  NPf'S  for  the  nex* 
iteration,  tut  ?OB?!C  aay  override,  e.g.  after  INvrPT. 

All  arrays  required  have  already  been  allocated  and  are  parasi'^ters  to 
DPIVER,  except  for  those  allccated  by  IhVEPT  for  the  Lb’  file,  naacly 
HI,  CL,  H'J,  DU  and  P. 

It  is  conceivable  that  CBIVSF  could  be  used  as  a subrouti..?  within 
sose  different  enviconaert,  hut  enormous  care  would  be  required  to 
ensure  that  the  various  arrays  and  labelled  COHHON  areas  were  set  up 
correctly,  in  the  way  that  is  presently  accomf  lisiied  by  IHITLZ,  SPFCS, 
SrPCS2  and  HPSIN. 


Parameters: 
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REALMS  Z(nAXZ)  (In, cut)  TV.®  array  of  coca.  IWTEPT  aast  b® 

ftfifl  to  usa  cor®  froa  Z (RZ1)  up  (t^aj^oraril y , 
'^hiie  factorizing  the  basis),  and  from  Z (KZl) 
up  thereafter  (for  the  LM  file),  where  KZ1,  FZ3 
are  preset  variables  in  CCHHOK  /COEF  /.  Of 
course,  all  the  array  parameters  cf  DFIVEP  also 
resida  in  Z in  the  current  i iplenent ation  of 
MINOS,  but  that  is  not  necessary  as  far  as 
CSIVFP  is  concerned,  and  need  not  po  »tii6  in 
scf3  cti.Qt  enviccnaent.  For  example,  tzi  and 
RZi  could  both  be  1. 


On  output,  the  solution  vector  begins  at 
Z(RZ5).  This  is  a vector  of  diaonsion  N,  and 
contains  values  for  strncturals,  [ 'S,  and 
logicals,  in  that  order. 

rWTEGBF  ZPSOLN  (Input)  Should  te  positive  If  ^he  solution  is 

tc  te  output  csing  the  standard  outpu*  routine, 
SCLS.  Should  be  zero  otherwise  (e.g.  if  the 
user  wishes  to  do  further  cowputation  on  the 
sc lut  icr.)  . 

IBTEGEBs  a,a2,HF,flNIi!,N,Ny,Jii''1,fJR,S»,NF0,nAXS,NS  are  all  input. 

Arrays  HE,RA,  A,  BI , EU,  HP , xij  are  all  input. 


Arrays 


HSPIKF,HS,PFlVL,'tPlVU,  FPIvr,BL,NtJ,B,C, 

X,PI,OBD,  T,  Y 1 ,G  1 , G.GciEW  pcovida  workspace. 


Paraaeters  «, H2, MS, RUN, N, HE, Sr 1 , MP, NN  ,NN0  ,MA XS  , 
HE, BA, A, Bl, EH 


will  be  unchanged. 


Host  other  paraaeters  say  be  used  as  output,  depending  on  the 
application. 


Method ; 

1.  Initialize  nuaerous  scalar  quantities.  This  includes  setting  the 
tolerances  XTCI(l),  FTGt(i),  gT0L(1)  to  their  "loose”  values,  as 
indicated  by  LVLTC1=»2. 

2.  Select  a starting  tasis,  with  the  following  order  of  pr^ifersnee: 
OID  BIT  HAP,  INSERT  JIl*,  LCAC  FILE,  CRASH. 

3.  Call  IMTERT  to  faetcrize  l.nltial  basis. 

4.  (Main  loop) 

(a)  Call  FORMC  to  test  teusjbillty  and  set  NPHS, 

(fc)  Call  CHKGRD  if  first  feasible  point  has  just  been  found. 

(c)  Call  BBSETR  If  the  Hessian  is  about  to  be  used  for  the  first 

ties. 

(d)  If  IPHS>2  and  feasible,  call  GB7GRC  to  get  an  appropriate  rhs 
for  coeputlnq  PI. 

(e)  Call  SBTPI,  unless  FT  has  already  been  coipetei. 
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(f)  If  IIPHS>2  and  HS>0,  call  CALCG  to  ccnpiite  the  redac»^d  qradient 


1 

L 


Jg)  If  llPHS<U^  call  PFICE  to  select  one  or  aore  additional 
siiperbasics. 

(h)  Skip  to  5 if  PBICE  found  nothing.  Otherwise,  reset  XT0L(1) 
etc.  to  theii;  "loose"  values,  .in  case  they  have  been  tightened 
up  by  a previous  skip  tc  5.  (Since  PtlCF  found  a variable  t« 
release  frca  its  found,  we  here  wrong  in  thinking  earlier  that 
the  final  set  of  active  constraints  had  been  detorsir  d. ) 

(i)  If  NFHS<=2,  call  LPITV. 

(j)  If  SPRS>=d,  call  EGITti.  On  exit,  if  NXTPHS  = 0,  PGITh  will  not 
have  perforaed  an  iteration  after  all,  but  wants  to  <-ry  again 
with  NPKS=1  (fcllcvinq  a PEICE) , 

(k)  Sisilarly,  a iteration  could  not  be  performed  if  IEPP=-1  (WTOPJ 
must  be  reduced)  or  if  INVPC>=20  (must  call  INVERT). 

(l)  Ctherstse,  an  iteration  has  been  successfully  cospleted.  Call 
ITEBOP  to  output  a line  of  Iteration  Log. 

(■)  Test  for  various  frequency  conditions  and/or  excessive  growth 
of  the  10  file.  Save  basis  and/cr  call  INVERT  if  necessary. 
Call  HOELU  it  reguested,  to  update  basis  factors.  Check 
residuals  every  KCFK  iteraticns. 

5,  (Apparently  optimal) 

PRICE  found  no  favorable  candidate.  If  infeasible,  sake  sure  that 

VTOBJ  is  zero;  If  feasible,  check  that  tolerances  XTOL(1),  etc.  are 

at  their  "tight"  values  (LVLTCL=3) . It  so,  we’re  probariy  done; 

check  residual?  or.3  last  tire  before  terminating.  Otherwise,  reset 

all  relevant  paranteters  and  repeat  from  u. 

6.  (Exit) 

(a)  Save  basis  map. 

(b)  Call  PDNGHD  one  last  time,  if  wanted. 

(c)  Set  up  solution  vector  in  Z(RZ1),  ...»  using  entry  1 of  SOLN. 

(d)  save  PUNCH  and/cr  DUHF  files. 

(e)  Output  solution  tc  printer  using  entry  2 of  SOLM  (unless  not 
wanted)  . 

(f)  Output  SOIUTICW  FILE  if  required. 
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SOBBOUTIPE  DOHEN(  H , H,  K£!  HS  , HS,  EL  , EO  , ¥ , IDl  , ID  2 ) 

Subroutines  called:  Called  by: 

Hone  CRIVEP 


Purpcse: 

To  output  a basis  description  to  file  IDOMP.  One  card  iaage 
containing  KEf,  HAHE  and  VALUE  is  output  for  each  variable,  except 
for  nonbasic  variables  that  are  fixed  cr  at  a zero  bound.  Logicals 
are  treated  the  sane  as  structurals.  It  is  intended  that  a DUKP  file 
be  easy  to  aodify  when  necessary  (that  is,  semewhat  easier  than  a 
PONCE  file) . 

Warning:  Nonbasic  variables  that  are  fixed  or  at  a zero  bound  will  not 
be  output.  If  the  file  is  used  as  a LOAD  file  for  a nodified  probles, 
it  nay  be  necessary  to  append  Li.,  OL  or  SB  cards  specifying  the-  state 
of  such  variables  explicitly,  if  their  hounds  are  relaxed.  Otherwise, 
the  initial  solution  obtained  aay  be  slightly  different  (and  in  some 
cases  fay  be  sericusly  infeasible).  Sivilar  croblems  arise  with  PUNCH 
fil  es,  but  they  can  he  overccsie  in  both  cases  with  a little  care. 


Faraeeters: 

PEAL*8  T (i)  (Input)  Contains  sclution  values  for  all 

variables  (structurals  and  logicals). 

THTEGBP  IC1(N)  Wor)c  vector  to  held  left  part  of  row  names. 

INTEGER  ID2(N)  Kcrk  vector  to  held  right  part  of  row  names. 


Method : 

1.  Bead  row  naees  fre*  the  scratch  file  into  IDI,  ID2. 

2.  Process  each  variable  in  turn.  Scae  eonkeying  around  is  necessary 
because  structurals  are  output  before  logicals.  For  each  variable 
d,  if  JCaKBHS  then  the  relevant  coluan  name  is  to  be  read  from  the 
scratch  file,  otherwise  tht  logical  nase  is  already  in  IDI,  ID2. 

1.  Skip  nonbasics  that  are  at  zero  bounds  or  are  fixed.  Otherwise,  use 
the  state  vector  HS  tc  index  the  local  array  KEY  to  get  the  correct 
indicator,  then  output  it  with  the  Naee  and  Value. 

u.  Print  sessage  to  indicate  successful  DUMP. 
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SnBPOOTlNE  FACTOF(  « , H 2, N , NE , NP  1, HF , « AXL A XO , 

1 HSP1KI,HA,BE,HE,HS,HPIVI PI» 0 ,H El V F , Hi , N 0, 

2 A,EI,Bn,HPIV9,HFIVI,X,Y,T1, 

3 Hl,DL,HU,D(l,F, 

4 TP1TR,TPIVC,T?'IN,HSKAE1,NS«AP2,NKELEM,KD0ME  ) 


Subroutines  called:  Called  by: 

BTFANL  FTHANL  PACKLO  UNPACK  INVEET 


Purpose: 

To  coapute  an  lU  f actorir ation  of  the  basis  aatrix  B,  as  specified  by 
the  list  of  coluan  nusters  in  HFIVO.  Arrays  HSPIKE,  HPIV3  and  HPIVI 
are  assuaed  to  be  set  cp  as  described  below. 


Paraaeters: 

I»TEGEP*2  HSPIKE  (H)  (Input)  Codes  tuap  and  spike  structure  as 

described  in  subroutine  P4. 

IKTEGEP*2  HPIVU  (n2)  (Input)  First  M components  are  ccluan  nos.  in 

the  preassiqned  order  from  P4, 

mTEGEP*2  HPIVP  (P)  (Input)  The  correspondinq  preassiqned  piwot 

rows. 

IHTEGEP*2  HPIVI  (H)  (Input)  The  inverse  of  HPIVfi,  used  to  test  if  a 

given  rew  lies  within  a given  buap. 

RIAL*3  X (M)  Work  vector  for  holding  a row  of  I inverse. 

REAl*8  Y (H)  Work  vector  for  bolding  transforiaed  columns. 

REAI*8  T1(K)  Work  vector  for  holding  pivot  elements  of 

alternative  spikes. 

BEAL*8  TPIVR  (Input)  Pow  tolerance.  Before  ccluan 

interchanges  arc  invoked,  a preassiqned  pivot 
element  is  tested  to  mate  sure  it  is  within 
TEIV"  of  the  bigoest  leaaining  urpivoted 
element  in  the  current  buap. 

REALtS  TPIVC  (Input)  Coluan  tolerance.  If  several 

spikes  are  available  for  coluan  swapping, 
the  one  finally  chosen  must  have  a pivot 
eleaent  within  I’EIVC  of  the  largest. 

PEAL*8  THIN  (Output)  The  saallest  row  pivot  ratio  that  was 

encountered  during  the  factorisation.  If  the 
factorization  has  to  be  done  again  for 

numerical  reasons,  TPIVF  aust  be  made  bigger 
than  THIN  ct  no  improveaent  will  occur. 

IWTBGEB  NSVAP1  (Output)  The  no.  of  times  a triangle  ccluan  was 

swapped  with  a spike. 

IBTEGEF  NSHAP2  (Output)  The  no.  of  tlaes  a spike  was  swapped 

with  scBC  other  spike. 

IHTEGE?  NBELEH  (Cutput)  The  no.  of  eleaents  of  L that  have 

been  imbedded  directly  in  A (those  in  non- 
spike coluans  if  IMBEP  = TPllE,  otherwise  0). 

IBTEGEP  KDONE  (Output)  A countar  to  save  werk  if  the  L'l 

factorization  has  to  be  redone  for  storage 

(rather  than  nuaerical)  reasons.  The  existing 
tow  and  coluan  ordering  will  be  retained  up  to 
the  pcint  KDCNE,  which  was  where  the  previous 
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factorization  got  interrupted 


Method : 


Gaossiaii  eliiaination,  vith  ccl'inn  interchanges  where  necessary  to 
preserve  nunarical  Interchanges  can  safely  be  restri'ted 

to  each  fcuap  separately.  Ir  necessaiy,  scbo  cclnans  Kay  be  re1act-»d 
froa  B and  replaced  by  si!ita.'i3  unit  vectors. 

For  K = 1,  ...»  H,  the  F-f.  cciuan  is  treated  according  to  the 
follcwing  steps, 

1.  Deteruine  if  K lies  vithi..  the  current  hunp  boundaries,  thus; 
IBGN  <=  K <=  lEKD.  I:  not,  use  the  bump  and  spike  information  in 
HSPIK?  tc  reset  IPG*J  and  :ENT). 

2.  The  K-th  basis  colasn  err  t ponds  to  varialile  J = HPIVU  (K)  , and  is 

suppoEOd  to  pivot  o'i  .O  i Tciv  = HPIVP  (K)  . Skip  the  column 
Iffaediately  if  .1  is  i s’ac’t  pivoting  cn  its  cwn  tow.  Otherwise, 

HSPIfiE(K)  is  negativs  if  che  crlumn  is  a spike. 

3.  (Son-spike)  The  coluirn  yill  bcccme  part  of  L directly,  unless  its 

pivot  elemaat  is  tcc  ssall.  A certain  "row  test"  must  now  be 
applied.  Ose  HPIV^  to  '"itAS,  the  largest  of  those  elements  in 

the  column  that  lie  in'?:,  +he  corrent  tump,  in  rows  that  have  not 
yet  been  pivoted  on.  ‘ .e  ccluisn  fails  the  tew  test  if  its  pivot 
element  is  saaller  ‘hat  IC.I  I?  'aa  absolute  test),  cr  TPIVH*TMAX  (a 
relative  test);  in  suci;  c.-cics,  skip  tc  5.  Otherwise,  copy  column 
J into  the  L-file  (if  I.'"'  DI  - cr  set  up  a pointer  to  the 

colusu  within  the  ccrs*-rai.i'  aatrit  A.  Then  repeat  ^rom  1. 

the  relative  ro-'  •<;.  f is  lot  applied  if  we  previously  had 
trouble  with  the  (-tb  crl'vs-!.  The  present  column  has  already  been 
determined  to  bo  the  , br'*itute  available,  and  GIVEOP  has  been 

set  to  TfOS  to  signi':/  •-ct. 

U.  (Spike)  ilnoack  ccluvr.  j n.  d transform  it  by  that  part  of  the  L file 
belonging  tc  the  currt'’i  lunp  (with  the  help  of  variables  LBGN  and 
JIBGN'  that  are  avalJahiv,.  to  FTFANl),  thus  obtaining  a vector  Y. 
Check  the  pivot  eleRen*-  is  » ?foi.s,  this  time  looking  at  elements  in 
T that  are  inside  the  currert  buup.  If  the  pivot  fails  tne  row 
test,  skip  to  5.  Cthervisr,  call  FACKLU  to  put  the  bottom  half  of 
Y into  L and  the  top  hi.lf:  ■■ '^to  n„  Repeat  from  1. 

5.  (Bow  tent  failed.  look  ''or  a rubstitwte  spike) 

(a)  If  GIVBOP  * TPl’F,  Ski;  to  • , «ie  have  already  tried  to  find  a 
suhstitu'' 0 for  this  r •..siticn. 

(b)  Scan  the  reaaining  sri-'cs  in  the  current  bump.  If  none,  skip 
to  7.  If  only  ciic,  grc.t  !''■  imaediatelv  and  skip  tc  6. 

(c)  Cthorwise,  find  the  It.TV-th  tow  of  (L  trans pose) -inverse , and 
use  it  ♦•o  ccppute  tt-  pijot  elenents  for  all  eligible  spikes, 
■^his  requires  eiErlo  Inner  rroducts  of  the  dense  vector  X with 
the  sparse  spike  cclvir.:.';.  At  the  same  time,  find  F"AX,  the 
largest  potentlnl  pirtt  ple?vnt  cn  row  IPIV. 

(d)  Now  go  through  the  ipile:  frcp  left  to  right,  skipping  those 


whose  pivot  eleiant  is  not  within  TPIVC  of  PRAX.  This  is  a 
"coluen  tast"  that  aakes  the  process  equivalent  to  Gaussian 
elialnation  with  ccluwc  interchanges.  If  the  current  column  is 
not  a spike,  take  the  first  spike  tha»  satisfies  the  ccluEn 
test,  because  it  is  likely  to  be  the  shortest,  and  skip  to  t. 
Ctherwise,  for  each  remaining  spike  that  satisfies  the  column 
test,  extract  its  preassigned  pivot  tow  (I.  PIV,  say),  and  find 
the  one  that  maxisizes  li  (IPIV)  . (This  heuristic  for  s&lectluq 
an  alternative*  spike  is  due  to  Charles  fCrabek  of  CDC.  It  is  an 
attempt  to  prevent  swapping  the  had  spike  into  a position  that 
will  again  require  a swap.) 

6.  Interchange  appropriate  elements  of  HPIVO  and  HSPIKE,  and  set 
GIVEIt  = TBUS.'  Fepeat  ftem  2, 

} 

7.  (Pe'Ject  coJ'ian)  Replace  cclutn  J by  the  ccrrespondicq  slack  (coluan 
no.  KPHStIPIV) . Accueulatfe  relacted  columns  in  array  JHEJ,  for 


SOBROOTINS  FORMC(  f* , N , NW  , N NO , N S , NS,  HB,  HS  , BL , BU,  X ,GiD,  XN  ,C,TOLX  ) 


Subroutines  called:  Called  by: 

FONGBD  CPIVEB 


Purpcse: 

i- 

1 To  deteraine  if  the  current  X is  feasible,  and  to  set  op  a gradient 

vector  G3D  to  be  used  ter  cenputing  PI. 

I 

I Also,  to  reset  the  phase  irdicator  NFH£  when  needed  (for  the  first 

^ iteration,  and  fcllcwir.g  IBVIRIS  . 


Paraaeters; 


REAl*8 


FEAL*3 


REALtfl 


REAI*B 


GBD (NS) 


XN  (NNO) 


C (NNO) 


TOLT 


(Cutput)  The  gradient  vector  for  the  current 
lineer  objective  (whether  feasible  or  not). 

cnDij)  = 0.0  if  X(j)  is  feasible; 

- 1.0  if  X(1)  is  above  its  upper  bnd; 

- “1.0  if  X (j)  Is  hulow  its  lower  bnd. 

GPD{ICBJ)  is  treated  specially  as  noted  below. 
(lr,out5  If  X is  feasible  but  the  previous 
iteratlcn  was  not  (or  if  INVEST  has  lust  been 
called)  , any  ncnbasic  nerlir.ear  variaoies  will 
Le  giv?c  the  appropriate  bound  values  in  XN. 
These  will  he  retained  during  subsequent 
iterations. 

(Catputj  Again,  or.  the  first  feasible  iteration 
or  after  INVERT,  the  gradient  vactor  for  the 
ronliteat  variables  Hill  be  ccwputed  in  C,  by 
ore  tall  to  FUNGPD. 

(larut)  The  feasibility  tolerance  for  all 
%atj.atles. 


Method : 

1.  Run  through  the  X(j)  vaijigs,  ccaputinq  the  nuaber  and  sua  of 
infeasibilities  (WIN?  and  3IHF)  and  setting  GPD(j)  as  above. 

2.  (Infeasible) 

(a)  If  this  is  the  fcsginning  of  the  tun  (ITN=0),  set  Nrns=l  as 

first  preference,  to  get  regular  Phase  1 siaplex  aethod.  Even  if 
there  are  sene  super t n.s  ics , NPHS  teaains  1 unless  NOLTIFLE  TRICE  is 
in  effect,  or  a basis  file  was  used  to  start  the  run.  (The  aim 
here  is  to  leave  superiasics  fixed  at  the  values  specifr‘=‘d  by  an 
INITIAL  bounds  set  auf"  pwrfota  necaai  Phase  1 - Phase  2 siaplex 

iterations  cn  the  reBair.xng  variables  as  long  as  possible.) 

(b)  At  this  stage,  GEC(IOed)  will  be  zero  unless  TARGBTina  is  in 
effect  (in  which  cass  X (ICBj)  aay  be  infeasible).  resetting 
GFD  (lOBJ)  to  GRD(IOBJ)  - N JNTES^VTOE J has  the  required  effect  for 
targeting  and/or  the  cooj.osite  objective  aethed,  and  also  gives 
zero  as  required  if  neither  is  In  effect. 


( 

f 
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3.  (Feasible) 

(a)  If  this  is  th«  fccgirniog  of  th*?  run,  set  mphs=2  for  cold  starts 
s where  possible,  for  the  reasons  given  in  2(a)  above. 

j ft)  In  all  Cases,  reset.  GFP(IOcJ)  = -PINIflz.  This  will  he  the  only 

nonzero  ccmponent  of  GFD,  and  gives  ths  correct  gradient  vector  for 
the  linear  objective.  G-TGBD  overwrites  any  nonlitear  can ponent*? 
later  on. 

(c)  If  this  is  the  first  feasible  iteration,  or  the  first  iteration 
after  IIVEBT,  set  up  ccnfcasic  nonlinear  values  in  XK  ari  coinput*' 
; the  nonlinear  objective  and  gradient  by  a call  to  FONGFD, 


I 


I 


^ 


m 


SOBBOOTINE  FTPKNK  M , H 2 , NE , NP  1 , P AXL  , HA  , H E,  Jl , H PIV L,  Ml, H L , DL,  T ) 
Subroutines  called:  Called  hy: 

None  ADCCOl  FACTO?  LPITK  PGITN 

SBTX 


Purpose : 

To  solve  the  systea  L*x  = y,  where  both  * and  y are  stored  in  the 
parasetet  Y,  This  Is  the  1st  part  of  a conventional  "FTRAN” 
operation. 

Faraaetecs: 

All  previously  defined. 


Method: 


The  column  ttansformaticns  constructed  ty  FACTOR  (during  INVERT)  are 
processed  first,  beginning  at  the  one  marted  by  LPGN  ar.d  JLDGN. 
During  INVERT,  t hes®  parameters  point  to  the  start  of  the  current 
bump,  since  any  earlier  transfer  rations  would  have  no  effect  on  the 
vector  being  transfer ved.  During  subsequent  iterations,  L3GN  and 
JLBGN  point  to  the  beginning  of  the  L file. 

If  NETAR  > 0,  some  additional  transformations  of  a different  •*-ype  have 
been  added  to  the  I,  file  during  iterations.  These  are  processed  next. 
Each  transformation  Is  really  a sequence  of  mere  el-‘mentary 
transformations  of  the  type 

"Add  DL  (J)  times  tow  I to  the  pivot  row  IPIV" 


where  I is  either  HI  (J) 
value  of  IPIV.  Thus,  a 
transf ermatnon  is  to 
stability  hy  the  method 


(if  that  quantity  is  positive)  or  the  previous 
negative  value  of  HL(J)  signals  that  the 
he  preceded  ty  an  interchange,  as  required  for 
of  Bartels  and  Golub. 


SUBPOOTINE  FTP*N0(  M , 1 2 , N F ,K A > U , HE I V U , H PIV? , KH , H D, DO, F , Y ) 


Subrcutines  called:  Called  by: 

Bone  CHOzg  ppics 


Furpcse : 


To  sclv3  the  systen  0 (transpose) ♦x  = y,  where  both  x and  y are 
stored  in  th®  paraireter  Y.  Ibis  is  the  1st  part  of  a conventional 
"PTFAN”  operation. 


f.ethcd : 


Porward-suhstitution.  Ir  this  case  0 and  V are  treated  separately. 
Coiutpiis  of  V (rows  of  U (transpose)  ) are  processed  from  fror.*-  to  hack. 
Ho  tcansforaations  can  be  s)cipped.  However,  since  P is  a d«nse  matrix 
it  is.  po.ssitle  to  traverse  either  by  rows  or  by  columns  as  convenient. 
Hers  the  forward-substitution  runs  across  the  rows  of  F and  thereby 
allows  a row  to  be  skipped  if  the  corresponding  component  of  the 
ccuputed  solution  is  negligible. 


i 

1 

■i 

I 

i 

i 
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SUBFCUTINE  "IINGfCi  N N , KS  , HE  , XN  , C , X,  F N ) 


1 

I 

.1 

Subroutines  called:  Called  by:  ^ 

1 CALCFG  AOCCCL  DPIVRP  POP^.C  FGITN  "1 

' SEAfiCI;  i 


Purpos-a : 

To  call  the  user-writton  subroutine  CALCFG  with  appropriate 
naraieters,  and  return  the  current  value  cf  th-^  obiective  r unction 
(both  linear  and  nonlinfar  terirs)  . will  I ? called  laly  if  the  current 
point  IS  feasible. 


Parameters: 


REAIV8 

XN(HN1 

(Input) 

Current  value  of  nonlinear  variables. 

a 

1 

} 

FEAL*8 

X (IS) 

(Input) 

Current  value  of  basics  and  super- 
basics. 

PEAL*8 

C (BN) 

(Cutput) 

Cemputad  gradient  c*'  nonlinear  oM. 

REAL *8 

FN 

(Cutput) 

Peturns  the  cbjectiva  value. 

1 

Fethod: 

1.  Extract 

basic 

and  superbasic 

values  free  X and  store  in  vk  in 

natural 

order. 

It  is  assuned 

that  nonbasic  nonl:near  variables  are 

already 

stored 

in  XN;  this  i 

s arranged  by  FOFFC  after  each  call  ♦•o 

1 N V E FT  . 

2.  Test  Mf'X  for  first  entry  to  FUKGFD;  set  parameter  ESTATE 
accordingly  and  call  CALCFG  to  obtain  nonlinear  objective  value  fn 
and  gradient  vector  C. 

3.  Test  aode  to  see  If  user's  subroutine  CALCFG  reguests  termination 
of  the  run, 

«.  If  a linear  obtective  rev  eiists,  change  fN  to  FB  - I{IOEJ)  to 
include  the  linear  part  of  the  objective. 
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SUBBCUTZHB  GZ1GFC(  MK  ,NNO  ,!•£ , F B,CrGFD  ) j 

Subcoutlnes  called:  Called  bf; 

lone  «OCCOL  DBIVEfi  RGITI  SIABCH 

Purpose: 

To  set  up  the  gradient  sactot  for  basic  and  superbasic  Tariables  in 
the  erder  defined  by  HE. 

Paraaetets: 

All  pi^eviously  defined. 


Method: 

GSD(i)  is  zero  for  linear  zatiables,  ♦C(k)  or  -c  (k)  for  nonlinears, 
where  k»HB(j).  GPD(ICBJ)  = -HZhiriZ. 


SUBROUTINE  GO  ( Z,N«CCBE  ) 


SubEcutines  called: 


M1N05 


Called  ty: 


MAIN  FROGFAH 


Pucpcse : 

To  call  MINOS  repeatedly  until  paraaeter  lEFRCR  is  negative,  which 
signals  that  no  furtheL  SPECS  exist  on  file  ISf^CS  and  hence  all 
protleas  have  been  pcocessed. 

This  rcutine  is  the  one  that  aay  be  altered  if  MINOS  is  to  b-a  used  as 
a subroutine  for  any  special  purpose.  For  ejcaarle  it  irav  include 
calls  to  a aatrix  generator  and  a report  writer  before  and  after  the 
call  to  KTHOS. 


Faraaeters: 


PEAL*8  2(NHC0PE) 

INTEGER  KWeeSE 


(Input)  The  available  array  of  core. 
(Input)  The  no.  of  words  of  core  in  7 


Method : 

The  standard  version  of  GC  is  shown  on  the  following  page. 
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SUBFOUTIIB  no  ( Z^NHCCBt  ) 
inFLICIT  BEAL*0  (C-G,0-Z) 

EFAL*8  Z(NilCOPE) 

CCMHON  /FILES  / I SCB,  INP  UT,  lOLDB  , INEW  B , IIISPT,  IPBCH,  ILO  AD  ,1DU(!P 
COMHON  /SCLMCB/  I SOIH, KSCLM , MSCl N , HST ATE, LCHKGR 


C SIAHDAHD  CALLING  ECUTINE  POP  »•  I N 0 £ 

C — 

C SOLfES  PBOPLENS  C«E  AL  A TI«E,  BEADING  SPECS  PBCH  CABD  F.IADE8. 

C STANCAFD  OUTPUT  BCOTINE  fEQUESrPD.  SHOWS, NCOLS,  ETC.  IGKORED. 

C 

c 

c 

100  ISPECS  = 5 
ISCFCH  « 8 

irsciN  * 1 

C 

C 

C 

C CALL  HA?GBN(  Z, N WCOPI , ISPICS , INPUT, lEFFOF  ) 

C IP  (lEPFOfi.NE.Cj  FEIUFN 

C 

c 

c 

CALI  HINOS  ( Z, NHCCFE, ISPECS, ISCFCH, IPSCLN, 


1 IEPnCF,NBCWS,NCOLS,lXS,LXL,LPI,LHS,LPREE,HFFFE  ) 

C 

c 

c 

C IF  (lERROF.NE.O)  GC  TO  9C0 

C CALL  REPBRK  ISCFCH,  ISCLN  ,NPOHS,  NCCLS, 

C 1 Z (LXS)  ,Z  (LPJ)  ,Z  (LFREE)  ,NFPBE  ) 

C 

C 

c 

900  IF  (lEFPOB.GE.O)  GC  TO  100 
BETOPN 

C END  OF  GO 

END 
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SUBBOOTINB  HASH  ( LF N ,N 2N , NCC ll , 

1 KBT1,KEY2,''OCE,KEtTAfa,NAMB1,  NAMR2,KA,FOUHD  ) 


Sobcoutioes  called:  Called  by: 

None  rFS 

Purpose: 

To  Icok  up  and/or  insert  entries  in  a table. 


Paranetors: 


INTEGER 

INTEGER 

INTEGER 

INTEGER 

INTEGER 

INTEGER 


INTEGER 

INTEGER 

INTEGER 

INTEGER 


LOGICAL 


lEN 

MEN 


NCOLL 
K2Y1 
KEY  2 
H03E 


KEYTAR  (LEN) 

HANE1  (NEN) 
NAf1E2  (NEN) 
KA 


FOUND 


(Input)  The  lor.gth  of  thd  has-h  table. 

(Input)  MdxiBiin  nuaber  of  entries  allowed. 

(lEK  >=  NFN.  Usually,  LEN  is  about  2*NEN.) 
Accumulates  total  no.  of  collisions. 

(Input)  left  pact  of  current  entry. 

(Irfut)  Right  part  of  current  entry. 

(Input)  1 for  Icot-up  only;  2 for  look-up  and 
entry  into  table  if  KTYI  and  KEY2  are  rot 
already  in  table. 

A table  of  keys  pcintim  in^o  the  list  cf 
distinct  ontries. 

Left  part  of  list  of  distinct  entries. 

Fight  part  of  list  cf  distinct  entries. 

(Output)  If  H0DE=1  (lock-up),  KA  points  to  the 
position  in  the  NANE1-NAWE2  list  wher-; 
KEY1-KEY2  was  found;  will  be  zero  if  r.o’" 
found.  If  (10D£=2,  KA  points  to  position  where 
new  antty  was  aade  (excep*  KA=0  if  table  is 
full)  Cl  to  position  wht re  existino  entry  was 
found. 

(Cutput)  TRUE  if  entry  KEY1-KEY2  already 
existed. 


Hethod : 

See  F.P.  Brent,  "Reducing  the  retrieval  tlae  of  scattf»r  storage 
techniques,"  Covi.  ACH  if>  ( 1973),  rp.  19t-109.  This  version  has  been 
siaplificJ  for  the  case  where  no  entries  will  ever  bri  deleted. 

Brent’s  aethod  is  well  suited  to  this  particular  application  because 
each  entry  (SCW  NAKES  in  the  «PS  constraint  data)  is  quite  likely  to 
be  lockfed  up  several  timet. 

Warning:  tha  hash  function  to  be  used  in  two  places  in  this  routine  is 
machine-dependent.  The  required  properties  for  it  are  docuaented  in 
the  dlBOS  User's  Guide,  srctioc  IX. 2. 
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SUBEOOTIME  IHITL7 


Suhtoutin«s  called: 


Koae 


Called  fcy: 


HIUOS 


Purpose: 

To  initialize  the  eaebine  precision  EPS  and  the  word-length  indicators 
NWOFCR,  NVOBDI,  NMOFDH,  and  to  conpute  nunerous  tolerances,  in  terras 
of  EES  where  applicable. 

Hethod: 

The  IBK  370  version  of  INZTLZ  is  shown  on  the  following  page. 


r>  o o o o o 


C 

C 

C 


C 

c 

c 

c 


SOBFOOTIIE  INI112 
IMPLICIT  EIAL*e  (C-G,0-?) 


CCHHON 

CCflHON 

COMMON 

COMMON 

COMMON 


/CJCCP.  / 
/EPSCCM/ 
/HGTCLS/ 
/TOLS  / 
/KOBDS2/ 


•:CLD.n,TCLC.I2,TOLDJ3,TOLDJ 
EPS,EF?0,EFS 1,PPS2,EPS3,FPS4,EPS5,PIINFY 
XIOL  (J)  ,rTCL  (3) ,CTOL  f 3)  ,P: NOPM, BG»OPM,TOLRG 
TCIX,TCLPIV,TEFI V1,TBPIV2,TOL?OH,XNOP« 
NWOPDF,NWORDI,NWCFDh 


THE  FOLLOWING  <4  NOMBIPS  APE  ThI  ONLY  MACHINE-DEPENDENT  FARAHETEB3 


c 

FPS 

THE 

MACHINE’S  FLOAT 

IKG-POIK 

T FFECISION 

C 

NWOBOB 

NO. 

Of 

"EEALS" 

PEB 

WOFD 

(VAF3  STABTING 

WITH 

A-B) 

c 

NWOBDI 

NO. 

OF 

"INTEGEFS" 

FEB 

WOFD 

(VAPS  STABTING 

WITH 

I-N) 

c 

NUOBDH 

9 

NO. 

CF 

"HALF  INTEGEES"  PEB 

WOPD 

(VARS  STABTING 

WITH 

H ) 

c 

HHEBE 

"wopr" 

BEANS  SPACE 

USED  BY 

VAPS 

STARTING  WITH  C 

-G,0- 

Z. 

IBM  360  AND  370 


FPS  » 16.0-**  (-13) 
NiORDB  « 2 
NWOBDI  = 2 
NHOHDH  = U 


USE  EPS  TO  SFT  CTFEB  HACBINE  PRECISION  CONSTANTS 


EPSC  = EPS** (4. C/5) 
EFS1  = BPS**(2.0/3) 
EPS2  = EPS**(1.C/2) 
EPS3  * BPS**(1.3/3) 
EFS4  = EPS**(1.0/4) 
EPS5  = BPS**(1.0/5) 
PLINFT  = t.0E»30 


SFT  TOLERANCES 


TCIX  * 

DMAX1  (EFS2, 

1.0C-5) 

TOLPIf  ■ 

FPS2 

TEPTVI  » 

0.00  1 

; TEPIV2  » 

0.1 

TOLBOW  ■ 

DMAX1  (EPS3, 

1.0C-4) 

; TOLDJ1  = 

DFAXI  (EES2, 

1 .OC-6) 

TOLDJ2  * 

1.0 

' TCLDJ3  = 

DMAX1 (EPS2, 

1.0D-6) 

XTOL(2)  = 

0.1 

XT0L(3)» 

TO  IX 

PTOL(2)* 

XTOL  (2)  *0.  1 

'i  FrOl(3)» 

XTGL  (3) **2 

GIOL(2>* 

DMAX 1 (Ef S4, 

1.0t-3) 

* GTOL(3)* 

DMAXI  (EFS2, 

I.Mp-7) 

I ■ C 

I PETDBN 

i C FKD  OF  INITIZ 
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SOBFOUTINE  INSFRI  ( H, K , NH , NNO , HN , M A X S , NS , HS , HB, B L, BO, X , XN , 1D1 , ID2) 


Subroutines  called:  Called  by: 

NKSRCH  E6IVEK 


Purpose: 

To  set  up  a starting  basis  using  data  on  file  INS8T.  This  routine  is 
irterded  to  provide  compatitil ity  with  ccEraetcial  systens  by  accepting 
data  in  the  foreat  ptcduceu  by  the  ccnventional  PONCH  comirand.  (It 
also  reads  a file  produced  by  subroutine  PUNCH  in  MINOS.)  This  data 
is  of  the  fora 

KSS  SAfEI  NAHE2  VALOF 

• here  KEY  aay  be  one  of  the  standard  set  IL,  IIL,  XL,  XU  or  else  SB 

to  indicate  superbasic.  VALUES  are  used  only  if  KEY  = SB. 


Para ae ter s 

INTEGER 

NS 

(Cutput)  The  final  no.  of  superbasics. 

INTBGEE*2 

FB (MM) 

(Output)  Mill  contain  a List  of  basic  and 
superbasic  ccluan  nos. 

HEAL*9 

X (MN) 

(Output)  XO)  will  contain  values  for 

...,  M+NS. 

INTEGER 

ID1  (N) 

Horksuace  to  hold  left  part  of  variable  naaes. 

INTEGER 

ID2  (N) 

Wcrkscace  to  hold  right  part  of  variable  naaes. 

flethod  : 

1.  Fead  row  and  coluan  naaes  free  scratch  file  into  ID1,  102. 

2.  Set  structurals  to  be  ncntasic  at  their  seallest  board. 

3.  Set  state  of  logicals  to  be  trasic. 

•.  Process  data  cards  cne  at  a tise.  Exchange  cards  are  ignored  if 
the  incoaing  variable  (NAME))  is  already  basic  or  superbasic,  or  if 
the  outgoing  variable  (KAKS2)  is  net  tasic.  This  guarantees  that 
there  will  always  be  H basic  variables  cn  exit. 


SOBEOOTINB  INVEST  ( « AXZ,  S , M 2 , M N , » , N t , NF 1 , NE  , MB  ,*J  ■ , «S  , K A X I , f* A XU  , 

1 HSFIKE, KA,FE,bI,HS,HPIVL,HFIVO,HFTVF, NL,NU, 

2 A,B1,BU,R, X,HEIVF,hFIVI, Y,Y 1,Z  ) 


Subicutiaes  called: 

DELCCL  FACTOR  F4 
SETX 


EEbETP 


Called  by: 


EFIVEB 


Purpose : 

To  cowpute  an  LIJ  f actcrizaticn  of  the  current  basis,  B, 


Paraioters : 

INTEGER^Z  I 
INTEGEB*2  1 
RFAL*fl  : 


PEAL+B 

BEAL*8 

REAL*« 


Method : 


HPIVR  ttcrkspace  to  hold  row  per3>iitation . 

HPIVI  (B)  Workspace  tc  hold  inverse  of  HPrVP. 

X (RN)  The  first  R locations  in  X are  available  for 

wcrkfpace  durinq  paCTOP,  but 
X (M  + NE)  contain  values  for  superbasics  and 
must  not  be  overwritten.  On  exit,  values  for 
basic  variables  will  be  in  X(1),  X(b). 

Y(H)  wcrkscace  for  FAr"'CP,  SETX. 

Y1(H)  Workspace  fcr  FACTCE,  SETX. 

Z(HAXZ)  Available  cote  beqins  at  Z(KZ1)  for  P4,  and  at 

ZIKZ3)  for  fACTOF.  (KX"!  < KZ3,  but  as  it 

happens,  FP  needs  iruch  less  core  than  FACTOP.) 


tse  the  state  vectcr  HS  to  select  column  numbers  of  basic 
vaciatiles.  Store  thase  backwards  into  PE,  so  that  slacks  will 
appear  first  and  qet  dealt  with  directly  by  the  transversal 
finder,  TES3VL.  Count  rerzeros,  slacks,  etc.  in  P and  print  one 
line  cf  basis  statistics. 

Call  P4  to  determine  the  buap  and  spike  structure  cf  B.  Host 
parameters  "just  provide  workspace.  Useful  output  is 
H£  iva  Bow  permtation 

HPIVn  Permuted  column  numbers 

HSPIKi’  Coded  bump  and  spike  structure. 

If  there  is  a linear  ob''€ctive,  make  sure  it  p'ivots  on  its  own  row 
'because  we  want  XdOEJ)  to  be  the  obfective  value). 

*f  IiPPEQ*!,  use  PSEIKE  to  display  th®  bump  and  spike  sTucture 
<■:;  ♦he  printer. 

♦ HPIVI  * inverse  of  HPIVR,  ready  for  FACTCP.  Initialize  pivot 
*ol«rances  TFIvp,  TFIVC. 

■•’.locate  core)  Allocate  storaqe  for  L and  U,  usinq  the 

• value  of  NSFIKF  and  the  values  of  KL  and  Ktl  for  the 

• V.  s LO  file  (which  usually  will  lavs  been  updatod  aaout  S'' 
r»  • r.  J Siiould  qiv*-  ^ fair  indicaticr  of  the  relative  sizes  of  L 


i ^ 


If  ? 
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and  n during  the  r^rt  SC  updates).  The  five  arrays  HL  and  DL,  Hd 
and  DO,  and  r,  will  occupy  the  fchole  o£  available  core,  Z (KZ3) , 
Z(HiXZ),  in  that  order. 

(a)  Storage  for  F is  allocated  first,  free  the  top  down.  The 
diaensioa  of  F is  taken  to  be  WAXSPK  ♦ N’EWSPK,  where  HAXSFK  is  the 
aaxiaum  no.  of  spikes  so  far  (in  case  of  aiultiple  a*t;jiiipts  to 
factorire  the  sa»e  iasis) , and  NEWSIK  is  estraated  tc  include 
about  10  percent  triangle  swaps  in  PACTCF  (each  ot  which  adds  cn.^ 
spike),  and/or  up  to  soaewhat  less  than  KTNV  nplates.  (No*  every 
update  results  in  an  increase  in  spikes,  since  a spike  say  l.e 
deleted  during  the  update.) 

(b)  The  resaining  core  is  allocated  tc  1 and  0 in  the  saj-^  ratio 
as  the  previous  LU  file,  lust  prior  tc  the  current  INVEPT,  naaely 
the  ratio  B = KU/KI.  Tci  safety,  this  ratio  is  meved  inside  the 
range  (0.2,  5).  At  the  ctart  of  a run,  hL  and  KD  are  initialized 
in  DFIVEF  to  previde  a cough  ci’ess  for  E as  follows.  s'or  a warm 
start  (anything  other  than  CFASH) , we  guess  that  L will  be 
somewhat  bigger  than  U and  therefore  take  KL=6,  KU=4  to  get  B = 
C.66.  For  CBASE,  although  L will  be  the  whole  of  the  initial 
(triangular  basis)  and  U will  be  I,  we  guess  that  the  rate  of 
growth  of  0 will  be  sutstantial,  and  therefore  set  K(J  = 2*FL  to 
gat  B = 2.0. 

7.  (Factorize  basis)  Call  FACTOR  with  the  current  list  of  basic 
columns  in  HPIvo.  (This  list  will  be  pertr.uted  by  FACTOP  if  column 
interchanges  occur,  and  will  be  altered  if  any  columns  are 
rejected  from  the  basis  and  replaced  ty  slac)es.)  Test  for  error 
condition,  which  can  only  cccur  if  insufficiert  storage  was 
allccated  for  L,  n or  F.  If  necef^sary,  repeat  from  6 to  try  a 
different  allccaticn  of  stcrag*:. 

(K.B.  The  Ionic  for  re-allccating  storage  for  L,  U and  F should  be 
iaptoved  for  the  case  where  tlere  ate  nary  spikes  and  not  much 
core.  At  present,  if  the  total  cote  available  is  only  marginally 
greater  than  absolutely  .necessary,  tbe  estimate's  may  fail  again. 
It  is  far  better  to  restart  the  tun  with  mote  core  for  Z.) 

8,  (Compute  basic  X)  Call  SETX  to  ccipute  the  basic  variables  X(1), 
...»  X (R) , using  ore  step  of  iterative  refinement  (without 
accumulation  of  residuals  in  extended  precision).  The  erj or  flan 


will  ho  St  if 

t he  row  check 

fails 

( i.  • Q • f 

if 

Ax  = b 

is  not 

sufficiently  well 

satisfied) . 

If  necessarv 

9 

t i ghte  n 

up  the 

tolerances  TPTVF, 

IPIVC  and  try 

agai  n 

f rria  b 

• 

TEIVF 

must  be 

made  smaller  than 

T.*i:n  (output 

£r  c B 

FACTOP) 

or 

the  pivot  order 

will  net  change. 

Initialize  arrav  RPivo  tc  be  HPIV' 

(they 

start 

to 

diffor 

during 

updates) . 

10.  If  any  variables  were  rejected  by  FACTO?,  check  for  suoerbasic 
slacks  that  ere  now  in  the  basis.  Call  DELC3L  if  necessary  to 
eliminate  them  ftem  the  superbasic  set. 
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SOBPOUTIRE  ITEFCF<  NS, NS  ) 


Sobroutines  calleil:  Called  by: 

None  TRIVEF 

Purpose: 

To  output  one  line  of  the  Iteration  log,  according  to  the  regnirfd  Log 

Frequency,  KLOG. 

Paraeeters: 

INTEGER  NN  (Input)  The  no.  of  ncnlinear  variabler. 

INTEGER  NS  (Input)  The  no.  of  superLasics. 

Method : 

1.  If  this  is  the  firrt  iteration  since  INVFBT,  initializ*.  COMFJN 
variable  IHi\D  tc  zerc.  (ItlAD  is  in  cCMKOh  so  that  it  will  retain 
its  value  between  entries.) 

2.  Exit  if  MOD(  ITS,  KLOG  ) is  nonzero.  Output  not  wanted. 

3.  Increeent  IREAD. 

4.  Branch  jf  proble®  is  nonlinear,  or  if  NPHS>2,  (Tore  inforwa* ion  is 
printed  per  line  fcr  non-sinplex  iterations.) 

5.  If  IhEAD=1,  this  Is  the  first  iteration  tc  be  prinned  sine-*  ■'hVERT. 
Print  relevant  Log  heading. 

6.  Print  one  line  of  the  Leg. 


SUBFOUTINE  LOACH  ( M , K , NH ,NNC , «N , HAXS , NS , HS , HB , BL , BO, X, XK  ) 


Subrcutla«s  called:  Called  by: 

None  CPIVEP 


Par  pose : 

To  input  a bit  aap  fro»  file  lOLCS  describing  the  state  ot  each 
variable,  along  with  the  cclunn  nuabers  and  values  for  superbasic 
variables.  (I.e.,  the  foraat  output  by  subroutine  SAVEP.) 


Parameters: 


INTEGER  NS 
INTEGEP*2  HS(N) 
INTEGEE*2  Hii(RN) 

BEA1*8  X(HN) 

REALV8  ZN(NSO) 


(Cutput)  The  no.  o£  superhasics  loaded. 

(Output)  The  state  vector. 

Supertasic  cols  HE(M^HS)  will  be 

o u t PU  t . 

Superbasic  values  X(M+1),  X (S+HS)  will  be 

out  put. 

Any  values  corresponding  to  nonlinear  variables 
will  t€  output  in  XN. 


Hethcd : 


Essentially  the  reverse  of  SAV5E,  except  that  several  consistency 
checks  aust  be  perfcroed. 

1.  Print  aessage  ‘BASIS  TC  EE  LOADED  FFOH  FILE  <IOLDB>'. 

2.  Bead  and  print  first  two  card  images.  Extract  M,  N and  NS  values 
fro*  the  second  card;  check  N and  N for  consistency  with  the 
current  probl-t®.  Ectct  exit  if  check  fails. 

3.  Read  the  state  vector  RS.  Exit  if  End  of  File. 


4.  Set  NS=0.  Load  coluan  nutbers  J and  values  XJ  for  superbasics,  one 

pair  per  card. 

(a)  Move  XJ  values  inside  bounds  if  necessary. 

(b)  Skip  if  HS(J)=3  (basic). 

(c)  Overwrite  HS (J)  with  the  value  2 (superbasic).  This  allows  the 
user  to  chanqe  ncntasics  to  soperbnsics  at  specified  valu<^s  by 
slaply  adding  cards  at  the  end  of  a sav<^d  bit  aap,  wi^-hout 
having  to  chanqe  the  PS  vector. 

(d)  If  the  value  XJ  is  essentially  on  a bound,  or  if  the  liirit  on 
superbasics  has  already  been  reached,  make  variable  J nor.basic 
at  the  appropriate  hound. 

(e)  Otherwise,  incienent  MS  and  store  J.  XJ  in  HB  (n«NS) , X(nvNS) 
respectively. 


I 


(f)  Stoc«  XJ  into  XN(J).  (Not  nocftssary  for  MINOS,  but  rtguirad 
for  BIMOS/GBG.) 


5.  Check  consistency  of  the  state  vector  HS. 


(a)  Check  that  nontasic  variahlas  are  not  specified  tc  be  at 
infiolte  bounds.  Switch  their  state  if  necessary. 


(b)  Count  the  nuaber  of  variables  of  each  state  2 and  1. 
exit  if  these  do  not  agree  with  MS  and  H. 


Error 
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SOBBOOTIWB  LOADS  ( M,  H,  KN  ,NliO  , MR  , MAXS  , NS  , HS  , HB,  Bl,  30,  X,  XN  , I r 1 , 1 1 ) 


Subroutiaes  called: 


NHSBCH 


Called  by: 


ERIVEB 


Purpose: 

To  input  a basis  description  froa  file  IIOAE,  in  the  fornat  output  by 
subrcutine  DUMP.  variaties  are  specified  by  state,  naae  and  value. 
This  type  of  basis  file  is  irtended  to  be  easier  to  modify  than  a 
FORCH/IRSERT  deck. 


Paraaeters: 

IKTEGEB  NS 

INTPGSE*2  HS|N) 
INTEGE3*2  HB(Mll) 

BSAI*8  X(NI) 


REALMS 

INTEGER 

INTEGEB 


Method: 


XN  (N«0) 

ID1  (H) 
ID2  (M) 


(Output)  The  no.  of  superbasics  loaded. 

(Output)  The  state  vector. 

Supertasic  cols  hB  (Me1\ , HB{WtNS)  will  be 

out  put. 

Supertasic  values  X(H*1),  ...,  X(aeNS)  will  be 
output.  The  first  M locations  of  X contain  the 
FHS,  cn  entry  and  cn  exit. 

Any  values  correspcnding  to  nonlinear  variables 
will  te  output  in  TN. 

Workspace  to  hold  the  left-  and  right-hand 
halves  of  the  row  and  cclunn  naaes. 


Straightforward,  except  that  aany  consistency  checks  must  be  applied 
to  ensure  that  exactly  M variables  end  up  basic,  that  no  infinite 
bounds  are  specified,  etc. 

1.  Print  Message  *LCAC  BASIS  BT  NAMES  — FILE  <IL0Ar>«. 

2.  Bead  and  print  the  first  card  iaage  froa  file  ILOAD.  This  should 
contain  the  pvotlea  caae  and  the  characters  ' DilMP/LOAD' . The  user 
can  check  visually  that  the  correct  file  was  loaded;  otherwise  it 
is  not  really  possible  to  deteraine  that  the  file  does  not  belong 
to  ECB€  other  problem. 

3.  Bead  row  and  column  naaes  frem  scratch  file  into  ID1  and  ID2.  Note 
that  row  naaes  are  first  on  th^  scratch  file,  hut  belong  at  the  end 
of  in  and  ID2. 


4.  Set  the  state  vector  BS  tc  aake  all  variables  nonbasic  ai 
siallest  bound  (in  absolute  value). 


their 


5.  Initialize  counters  and  proceed  to  read  and  process  cards  ori.'>  at  a 
tine,  until  an  ERTATA  card  is  found,  cr  End  of  File.  Each  card 
contains  a KEY,  NAME  and  VALUE. 

(a)  Call  NMSKCH  to  deteraine  a coiuan  no.  J corresponding  to  NAME. 
Skip  if  not  found  ( NMSFtH  will  have  printed  an  error  message 
and  incremented  the  error  counter  ). 
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(b)  If  J Is  a slack,  convert  cc«  value  XJ  to  logical  value,  using 
the  RBS  value  stored  in  X.  If  necessary,  sove  XJ  values 
inside  bounds. 

(c)  If  KEY  » 'ES*,  sake  the  variable  basic.  Skip  if  it  already  has 
state  3,  or  if  M basics  have  already  leen  specified. 

(d)  If  KEY  = ’LL'  or  'OL',  sake  tha  variable  nonbasic  at  the 
appropriate  bound.  Switch  state  if  necessary  to  avoid  infinite 
bounds. 

(e)  If  KEY  * 'SB*,  lake  the  variable  superbasic  as  long  as  it  is 
not  already  at  state  2 or  3.  Skip  if  the  no.  of  superbasics 
has  already  reached  MAYS  (tha  variable  will  reaaia  nonbasic). 

Print  no.  of  basics  and  superbasics  specified,  the  no.  cf  cards 
read  and  the  no.  ignored. 

(Partial  basis)  If  fewer  than  N basics  were  specified,  proceed  to 
wake  loglcals  basic,  startiro  frci  the  left  and  skipping  any  that 
are  already  basic  or  supertasic.  This  will  not  necessarily  give  a 
good  final  basis.  Actually  it  will  fail  if  toe  many  logicals  were 
already  saperbasic  — there  will  not  be  a full  set  of  r variables 
at  state  3 on  exit.  The  only  safe  way  would  be  to  add  a loop  to 
run  through  the  superbasic  list  backwards,  changing  as  many  states 
tc  3 as  necessary. 

Check  that  the  logical  cn  the  objective  row  is  basic.  If  not,  swap 
it  with  the  last  basic  variable. 


SUFFOnTIIIE  LPITN(  B , N 1 ,n  2 , N , HE,  M F 1 , N F,  P A XL  , HA  XU  , NM  , 

1 HA,HB,HE,HS,HPlVl,HPIVU,HFI7F,JiL,Nl),  A,BL,3(J,PI,1,T,Y1 

2 HL,DL,HO,DU,f  ) 


Sabroatines  called:  Called  by: 

BTPAHU  CHHZB  FTBANL  UKFACK  CFIVEE 


Purpese: 

To  perfora  a nornal  siaplex  iteration  after  PRICE  has 
variable  JQ  to  eater  the  basis. 


Paraaeters: 


INTEGER 

H 1 

PEAL*8 

PI(H) 

REAL*8 

X (HI) 

INTEGEB*2 

HB(H1) 

REAL*8 

Y (H1) 

REALtS 

Y1  (M) 

(Input)  Mel. 

Kot  used. 

(In, out)  The  first  tl  locations  contain 
for  the  basic  variables.  X(ni)  will 
for  the  iccciinq  variable. 

(In, out)  The  list  of  basic  variables, 
will  sieilarly  be  used  for  the  inccainq 
fccrkspace  to  hold  the  search  direction 
will  be  the  updated  coluan  A (JQ) ) . 
(Cutput)  will  contain  the  partially 
cclufltn  S (JQ)  . 


Method : 

1.  Set  logical  ATBND  tc  indicate  whether  variable  JQ  is  coaitiq 
its  upper  bound  or  not.  The  sign  of  DJQ  tells  the  story. 


se  lected 


values 

be  used 

HB  (Ml) 
var. 
(which 

updated 


in  froa 


2.  Set  XJQ  to  the  appropriate  bound  value.  Must  be  zero  if  JQ  is  a 
free  variable. 


з.  Unpack  colaar.  A (JQ)  intc  Y.  Solve  B*Y  = A (JQ)  by  solving  L*Y1  * 
It  U*Y  = Tl.  The  interaediate  vector  Y1  is  required  on  exit  by 
MODLU  to  update  the  LO  factors  of  3. 

и.  If  infeasible  and  a conposite  objective  is  being  used,  aodify  DJQ 
to  be  the  reduced  gradient  associated  with  the  su»  of 
infeasibilities  (uncentasiratf d by  a aultiple  W70BJ  of  the  real 
objective).  Exit  with  error  flag  set  if  the  aodified  DJQ  indicates 
that  the  sua  of  in f ea si tilities  would  net  decrease.  (DFIVEF  will 
reduce  WTOBJ  and  try  again.) 


5.  txit  with  an  IHTEPT  SEQUE3T  if  the  LU  file  hae  essantlally  no  space 

left  for  CHOZF.  (Actually,  CBOZR  does  not  need  space  if  the 

sclution  IS  feasible,  but  INVERT  would  be  needed  pret-y  scon 
anyway.)  Otherwise,  save  the  current  values  of  HF(Pi)  and  X(H1), 
which  aay  belong  to  scae  genuine  superbaslc,  and  install  the 
corresponding  values  for  variable  JQ. 

6.  Call  CHUZF  to  select  JP  -•  the  variable  leaving  the  basis  will  be 

the  JF-th  variable  in  the  list  HB.  Exit  If  JP-0  or  THETA  is  very 
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large  (problem  is  uotcunded)  . Treat  JF-*i  1 specially  — this  means 
the  inccilng  variable  has  reacted  its  of posite  bound.  Soaetiaes, 
ChUZe  mill  return  JP<0.  This  aeaas  that  the  (-JP)~th  variable  is 
currently  infeasible,  and  vill  become  feasible  this  iteration  when 
it  leaves  the  basis.  Record  this  fact  and  reverse  the  sign  of  JP. 

7.  Perform  housekeeping  to  update  the  list  of  basic  variables  and  to 
reset  the  state  vector  HS.  Heccapute  THETA  exactly,  in  case  CHUZR 
used  a perturbation, 

8.  If  TiETA  > EPSO,  update  X to  become  X ♦ (or  -)  THFTA^T.  The  value 
fcr  the  inccming  variable  is  then  in  X(H1|. 

9.  Restore  the  original  values  of  h.  (Bl) , Xtdl)  and  exit. 
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SnBFOOTIIIE  MINOS  ( Z , NWC06 K,  ISPECS , I SCPCH  , IPSOLN  , 

1 :FN»OF,NEOHS,NCCLS,[.y3,IXL,LPt,L«S,LFFEr;,NFfPB  ) 


Subroutines  called:  Callr:!d  by: 

OBI?ER  INItlZ  HPSIN  SPECS  GO 

SPECS  2 


Pucpcse: 

To  request  input  of  the  SPECS  file  and  MPS  file,  and  to  call  the 
solving  routine  DBIVEF.  Ihe  argunert  lists  in  the  calls  to  ether 
routines  here  Icok  pretty  ugly,  but  the  actual  routines  baing  railed 
(in  particular,  HPSTN  and  PPr?E?)  are  reasenabiy  ncrsal.  F*  wever, 
groat  care  is  required  if  ever  anything  is  aedified. 


Paraae'^ers: 


?EA1<8 

TMTEGEB 

INTEGER 

INTEGER 


IMTEGEP 


INTEGER 

INTEGER 

INTEGER 


INTEGER 


INTEGER 

INTEGER 


Z (NNCCRB) 

NNCORE 
I SPEC? 

IPSOI.N 


IRBROR 


(Output)  The  array  of  core,  us«d  for  all 
vcrlcEpace  and  to  return  certain  array 
infer  sation. 

(Inpiut)  Distension  of  7. 

(Input)  Onit  no.  from  which  the  SPECS  file  is 
to  te  read. 

(Input)  Should  he  positive  if  the  solution  is 
to  he  output  using  the  standard  output  rcunine, 
SCLN.  Should  he  zero  otherwise  (e.g.  if  the 
user  wishes  to  do  further  computation  on  the 
solution) . 

(Output) 

-1  if  FOE  occurred  while  trying  to  read  from 
file  ISFEC3  (signals  end  of  problems). 

0 If  optimal  solution  was  found. 

1 if  ptoblen  was  inteasitle. 

2 if  problen  was  unbounded. 

3 if  iteration  limit  was  exceed ^d. 

4 if  iterations  were  terminated  by  seme  other 
error  condition. 

30  if  there  was  net  enough  core  to  inpu*  the 
MPS  file  (NBCOFB  too  small)  . 

40  if  seme  other  fatal  error  occurred  during 
input  of  the  HFS  file. 


NBONS 

NCOIS 

LXS 


IXL 


IPl 

LHS 


(Output)  No.  of  roms  in  the  constraint  ma*:rix. 
(Cutput)  No.  of  cols  in  the  constraint  matrix. 
(Cutout)  AddC(=s8  in  z cf  the  sclutlon  vector. 
The  structnrals  form  a vector  of  length  NCOLS, 
starting  at  Z (LXS) . 

(Output)  Address  in  7 of  the  slack  variables. 
These  form  a vector  of  length  ■■fROfS,  smarting 
at  Z(LXL).  (Alternatively,  thw  strurtu~als, 
FHS  and  slacks  form  a vector  of  length 
NCC1S*1*NR0WS,  Starting  at  Z (LXS) . ) 

(Output)  Addtsss  In  Z of  the  dual  solution 
vector,  PI. 

(Cutput)  Address  in  Z of  the  sta**?  vector  HS. 


Hu 
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IHTEGEB 

IKTEGEP 


LfBEE 

KFPES 


This  is  an  IKTEGEP*;  uriny  of  ;JCaLS+l» 

KPCIi3  . 

(Output)  Eldcess  in  7 r,f  th*  'iiMt  iTftf* 
Iccation  in  Z- 

(Output)  Ko.  of  £re«  in  7..  Words 

Z(LfKEE),  ...,  Z(HUCOFE)  fliq  not  U8«d  to 
return  any  of  the  at.ovo  solution  vpttorR,  and 
h«jDce  nay  be  used  as  workflpi*c«  by  th**  ral  lin<j 
routine. 

If  necessary,  various  ether  address  pointers  couid  he  addad  os  output, 
parameters.  Those  pcioting  to  the  constraint  matris  and  bounds  are 
the  most  likely,  namely  rii £,  KH  * , RAX ,llf Bl , K RU . 

Eethod : 

1.  Call  INITLZ  to  initialize  a few  CCMMON  variables. 

2.  Call  SPECS  to  look  for  a SPECS  file  and  output  it  in  nperi.tl  fornat 
tc  the  scratch  file,  if  no  SPECS  file,  exit  with  lEPPOP  -1. 

3.  Call  SPFCS2  to  input  the  SPECS  paraaeters  froti  :icritrh  file. 

This  returns  various  diBiecslon  parameters. 

0.  Call  MPSIf  to  input  the  ftPS  file  and  to  allocate  u’l  array  storage. 
Pest  of  the  parameters  return  starting  addtesso.ii  within  Z,  for  the 
arrays  used  by  DEIVEB. 


SDBBOOTINE  KKIISl  ( W , H 1 , NB,  UP  1 , N Z,  H A , ti£ , HB  , li' , IP M ) 


Subtcutlnes  called: 
None 


Called  ty: 
FU 


Purpose: 


To  set  up  a colacm  list,  containing  the  positions  of  nonaeros  in  the 
set  cf  coluans  defined  ly  artay  HB. 


Parameters: 


INTEGER 

INTEGER 


(I 

HZ 


IKTEGRB*2  1P(H1) 
lNTEGEa*2  IFH  (MZ) 


(Input)  E*1. 

(Input)  Total  no.  of  nonzeros  in  the  list. 
(Output)  IP(j)  eill  point  to  the  start  of 
the  j-th  column,  in  IPN. 

(Output)  The  cciumn  list. 


Method : 
Obvious. 


S03FO0TINE  NODLO  | H , B 2 , VF , HAXl , KAIO , HPI V L, HPIfU, HPIf F, NL , MO , 
1 BL,CL,«U,rU,P,HE.X,Y,Yl,JP  ) 


Subrout-inea  called: 


Purpose; 


Called  by: 


caivEP 


To  Bodify  the  10  factctl7aticn  of  the  tasis  vrhan  the  JP' th  cclumn  is 
replaced  fcy  soae  vector  S (JC)  » usieg  Gaussian  elimination  with  row 
interchanges  to  maintain  numerical  stability,  a ia  Bartels  and  Gciub. 


Paraaeters: 


IN?EGFF*2  HB(K) 


RBAL*8 

HEAL*8 

REAL*8 


Bethod : 


Y(H) 
Y1  (H) 


(In, out)  List  of  basic  variables.  *’ust  be 
permuted  to  match  any  row  interchanges. 

(In, out)  Valv,js  of  basic  variables.  Same 
ccmnient, 

hcrk  vector  to  hold  the  JP-th  row  of  0. 

(Input)  Contains  the  solution  of  the  system 
L^YI  - A (JQ)  . Gets  overwritten. 


That  of  Battels  and  Golub,  as  implemanted  by  Saunders  (1976) . The 

tolerance  TOLSMP  (O.C  <=  TOLSWP  <=  1.0)  is  used  in  the  test  for  row 
interchanges.  1.0  gives  the  original  Bartels  and  Golub: 
simulates  Forrest  and  Tcmlia.  In  practice,  0.99  or  0.9  is  pl<jnty  big 
enough  to  preserve  stability,  and  reduces  worK  very  slightly  by 
avoiding  a few  interchanges. 

1.  Increase  HSPIRE  and  make  JP  a new  spike  row. 

2.  Use  HPiyp  to  pick  oat  all  spike-row  elements  from  T1  and  for*  a 
new  coluiin  of  F.  Reset  such  elements  of  Yl  to  zero,  so  they  will 
cot  be  packed  up  into  0.  Save  pivot  element,  which  will  become 
part  of  the  JF-th  tew  of  [ U Y 1 ]. 

3.  Pack  remaining  nenzeros  in  Yl  to  form  a new  column  of  U. 

h.  Test  if  the  JP-th  column  of  B is  a spike.  If  so,  the  JF-th  row 
of  "U"  comes  feem  F,  otherwise  from  U.  The  <:e8t  involves 
jumping  out  of  a loop  on  JO.  if  liPIVF(JC)  = JP-  This  test  will 
always  be  satisfied  at  least  the  last  time  round,  for  J0=NSPIKR, 
because  of  1 above.  Hence,  progcaa  will  not  fall  through  the 
loop. 

5.  Sat  Y equal  to  the  JP-th  row  of  "0".  If  a spike  is  not  being 
deleted,  this  means  runring  through  HU,  co  looking  for  ronzeros 
in  tow  JP,  which  ate  physically  reset  to  zero.  Otherwise,  JQ 
points  to  the  correct  row  of  F.  Install  this  row  into  v while 
deleting  the  J0~th  row  and  column  from  F.  Also,  JU  pcin*s  to 
the  spike  being  deleted;  make  HDivn(Jl))  negative  to  ii  M -ate  to 
FTSAHU  and  BTRANU  that  this  column  has  been  deleted. 
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6.  At  this  stage,  F is  alaost  upper  triangular,  except  for  its 
tottoe  row,  which  is  stored  in  Y.  Now  eliuinate  any  norzeios  in 
Y.  This  will  be  straightforward  Gaussian  eliaination  with  row 
interchanges.  The  pivot  rcw  for  Y <in  IPVIY)  is  initially  JP. 
All  pivot  rows  (IPIVY,  ariVF,  HFIVD)  will  reeain  unaltered 
unless  a cow  interchange  is  required. 

7.  Scan  Y forwards  locking  for  a nonzero.  Suppose  Y (j)  is  found 
to  be  significant.  Ccmpare  with  the  correspondir.g  diagonal  of  P 
as  follows: 

if  abs(  diagonal  cf  F ) < T01SWP*abs(  Y(j)  ) 

then  interchange. 

8.  If  no  interchange,  add  appropriate  wultiple  of  row  of  F to  Y. 

9.  Cthccwise,  interchange  Y with  the  row  cf  F and  perform  the 
elifination  simultaneously.  In  this  case,  the  pivot  row  for  F 
<and  0)  is  swapped  with  IFJVY.  Also,  the  relevant  ccmponents 
cf  HP  and  X ace  .i  ntercbangsd. 

10.  Pack  the  multiplier  and  pivot  row  (negative  if  an  interchange 
cccurred)  into  the  I file,  and  repeat  from  7 for  the  next  nonzero 
in  Y. 

11.  If  any  elements  of  Y were  eliminated  (often  there  will  be  none), 

increment  NFTAB  and  fill  in  NL,  FPIVL  for  the  new 

t rans forma tic n. 

12.  Test  the  final  diagonal  cf  F.  Request  INVERT  if  smaller  than 
PPSI  (IPS**  (2/1) ) . This  test  should  really  be  relative  to  the 
norm  of  A(JC),  but  we  assume  that  all  columns  have  approximately 
unit  norm. 
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SUBBQOTIIE  nOVE(  FI, R2,K«N 1,M2  > 


Sabroutln«s  called: 


Called  bj: 


HPSIN 


Purpcsa: 

To  copy  array  B2  into  a specified  pact  cf  array  F1.  This  is  required 
by  SPSIN  during  input,  once  the  diaenctcns  of  the  probletp  are  known. 
In  particular,  the  bounds  foe  the  slacks  have  to  be  ipoved  into  theit 
final  position.  To  do  this,  a subroutine  call  is  required  in  order  tc 
maintain  alignment  of  EfAl  bounds  within  the  FEAL*8  array  Z that  HPSIN 
has  to  deal  with. 


Parameters: 

BEAL  1 

BEAL  1 

INTECES  I 


net hod: 


B1  (Ml) 

F2(M2> 

B 


(Output)  B1(K*1),  ...,  H(K  + N2)  mill  contain  the 
contents  of  F2. 

(Input) 

(Inpat)  The  required  offset  in  Bl. 


Trivial,  as  follows: 

CO  10  J=1,?l2 
10  P1(K»J)  » B2(J) 
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SUBBOOTIHE  BPS  ( KCMl, 

1 LP«,NEN  ,IICCLI,ME0«S,BC01S,  MELBS,  *,  HE  , MP  1 , MB  , MNO  , NS  , *5  MS  , 

2 KCAFD,HfTYPE,kAfi?1F,«»MEia,l<ErHAR, 

3 BE,HA,A,SE,BO,HE,  Xii,N»R£1C,NABl2C  ) 


Sabroutinas  called:  Called  by: 

HASH  NMSRCH  BP  SIN 


Purpose: 

To  input  constraint  data  in  standard  BPS  foraat.  This  is  the  foraat 
used  for  specifying  large-scale  linear  prograBiing  proMeas  (the 
"COHVEFT"  data  formal  cf  T^f’s  HPS/360  and  Brsx/37C)  . 

nonzero  elaaents  in  the  corstcaint  oatrix  A are  specified 
colupn-aise.  Ail  elcacr.^s  in  a cC'Urn  lust  te  specified  on  ccnsectiTe 
cards  (there  is  no  ettcr  check  to  detect  split  ccluens) . 

Three  entries  are  required  tc  input  various  sections  of  the  data,  as 
specified  by  the  paraaeter  NCAll,  thus: 


NCALl 

* Incut 

NABE  and  PCWS 

2 Input 

CClilMNS,  PHS  and 

BADiGES 

3 Input 

EO'JNCS  and  IMTIAL 

BOUNDS 

The  FAKGES 

and  BOONES  sections 

need  not  be  present. 

Data  should  be 

terainated 

by  an  EKTATA  card. 

Paraaeters 

INTEGER 

NCALL 

(Input) 

Defined  above. 

INTEGER 

lEN 

(’'nput) 

Length  cf  hash  table  for  tov  nares. 

INTEGER 

NEN 

(Incut) 

No.  cf  distinct  entries  in  hash  table. 

INTEGER 

SCOLL 

(Input) 

No.  of  collisions  during  searching  of 

hash  table. 

INTEGER 

BBOtiS 

(Input) 

BaxiauB  no.  of  rovs  allowed  for. 

INTEGER 

BCCIS 

(Input) 

Baxifua  no.  cf  cols  allowed  for. 

INTEGER 

BBLHS 

(Input) 

HaxiauB  no.  cf  aatriz  eleaents  allowed 

for  . 

INTEGER 

W 

(Cut, in) 

The  actual  no.  of  rows  found. 

I>JTEGEP 

NCAEO  (fc) 

(Output) 

No.  of  cards  in  each  section  cf  data. 

INTEGER* 2 

HPTYPE  (BECKS) 

(Cut , in) 

Codes  constraint  types  as  follows: 

-1 

G >= 

0 

E 

1 

1 <» 

2 

N Non-binding 

IHTEG»P*2 

HAHE1B  (BECKS) , BABE 28 (BPCMS) 

(Cut  ,in) 

left  and  right  halves  of  row  r.anes. 

INTEGER*2 

HABEIC (BCCI 

S)  , NA''r2C  (BCCIS) 

(Gut, in) 

Left  and  right  halves  of  col  naaes. 

INTEGER 

KBTKAB  (LEN) 

Cut,  in) 

The  hash  table  pointers. 

PEAL*3 

XI  (NNO) 

(Output) 

After  entry  3,  returns  initial  values 

for  rcr.linear  variables. 
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Ret hod : 

The  required  OBJ,  RHS,  CANGE  aod  BCOID  Bases  are  in  the  /RPSCOR/ 
arrays  ROBJ  (2) , MPEG  (2),  RENC{2).  If  these  naass  are  blank, 
the  first  naaes  enco>jntece<^  in  the  rfs  data  will  be  used. 


Several  other  paraeaters  reside  in  labelled  COHRON  areas.  In 
particular,  the  variables  in  /RPELCl/  are  net  used  by  any  other 
subroutine,  but  are  in  a CCRRON  area  to  ensure  that  they  retain  their 
values  between  entries  to  RES. 


HCAll  * 1 


Input  MANE  and  PCMS. 


1.  Initialise  various  counters.  Zero  the  hash  table  pointers. 

2.  (hAHE) 

(a)  Fead  first  card  froa  file  INPOT.  Error  aessage  if  it  is  not 
the  RARE  card.  Extract  naac  froa  it  anyway  and  proceed. 

(b)  Skip  if  first  card  was  the  PCHS  card,  otherwise,  read  second 
card.  Error  aessage  if  not  the  ROiiS  card.  Proceed  as  if  it 
were  anyway. 

3.  (COWS) 

Read  cow  types  and  naees  until  the  CCLORNS  card  is  found  (cr  EOF). 

(a)  Allow  types  1,  G,  E,  M tc  cccor  in  colunn  2 or  colunu  3. 

(b)  Hake  IC3J  point  to  the  appropriate  N row. 

(c)  Enter  each  new  row  naee  into  the  hash  table  and  check  for 
dupl Icates. 

(d)  Hrite  row  naaes  cut  to  the  scratch  file. 


NCALL  = 2 


Input  CCIOHNS,  BH$  and  RANGES. 


1.  (COLOHNS) 

Read  data  until  RRS,  BCUNOS  or  BMCATA  cards  are  found.  For  each 

card: 

(a)  If  the  coluan  nane  is  different  froa  that  on  the  previous  card, 

increaent  the  coluan  counter  N,  save  the  new  naae  in  core  and 

also  write  it  cut  to  the  scratch  file.  Error  aessage  if  the 
cid  coiuwn  had  no  eienents,  unless  it  belongs  to  a nonlinear 
variable;  for  ncniineats,  enter  a duany  zero  in  row  1,  so  that 
Fortran  loops  car  go  round  once  without  damage. 

(b)  Extract  0,  1 or  2 rew  naaes  and  aateix  coefficients.  Lock  up 

any  non-black  rew  naae  in  the  hash  table.  Error  aessage  if  not 

found. 

(c)  Count  up  total  eleaents  EE,  Ignored  eleaents  RAO,  etc. 
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. (FH3) 

Very  sisllar  to  CCLUPI'S  irput,  except,  for  soae  iugqllng  tc  reiect 
cne  cf  poasibly  lultir^e  fHS’r.  If  no  PUS  eleaents  are  fcurd,  or 
if  all  nre  aero,  insert  a dMaiBy  FHS  vith  zero  in  ro*  1.  (This 
serees  ■»«  a place  hcldec  for  roluin  no.  KFHS.) 

3.  (PUNGFS) 

(a)  Set  defa'ilt  toundu  on  tha  slacks,  accordina  to  the  row  types 
(stotad  in  hrtyF!  ilurino  entry  1). 

(t)  Eilt  if  there  is  ac  r.ir.ge  set. 

fc|  Paad  cards  nntil  required  FANGB  riaie  Is  found.  Extract  row 
n»aef  and  range  va’ues  in  similar  fashion  to  CUI.(Jf!W5  and  PHS 
1 nr  I. » . 

M)  Treat  raroe  waloes  according  to  the  tow  types. 


aCALI  1 Input  PCOKC?  and  INITIAL  BCriNOS. 


1.  Aid  i fall  j.daatity  tattix  to  the  and  of  the  constraint  aatrix  (one 
ccrzerc  per  black). 

2.  Set  bounds  on  structnrals  tc  their  default  values  ESTPUC{H  and 
BSTPOC (2) . Thase  are  usually  C.C  and  1.0Et3C. 

3.  Skip  if  the  FNCAIA  card  has  already  been  read.  Ctherwise,  read 
data  until  the  required  bounds  set  is  found,  or  the  reserved  naae 
INIIIAL  is  encountered. 

4.  ’'cr  notaal  bounds,  read  one  bound  typo,  ccluan  naae  and  value  per 
card.  Do  a linear  search  to  find  which  coluar.  it  is,  starting  fro« 
the  point  where  the  previous  search  finished.  (It.  is  likalv  that 
bcurds  will  be  roughly  in  natural  order,  so  the  linear  search  w'll 
usually  teraiuate  very  quickly. 1 

5.  Pofoce  reading  any  l.MTIAL  bounds,  set  XN  to  the  ssallest  bound  (in 
dtsclite  value)  cn  the  corresponding  nonlinear  variables.  Then 
read  through  ail  ceaaining  data  cards,  until  EhtiTA  is  founi.  Evon 
if  the  bounds  are  to  be  ignored,  the  cards  j'ust  be  read  in  case 
further  problets  occur  later  in  the  input  file. 

6.  Court  the  no.  of  cows  and  cclusns  of  various  type  (MOBHAL,  FFEB, 
FITEC,  BOUEPBD  and  output  as  flAlFIX  STAIISTTCS.  In  this  centext, 
NCF^AL  aeans  variables  with  bounds  0 to  INFINITY  or  -INFINITY  to 

0. 


97 


SUSaoatlNE  «PSIN(  Z,H«XZ,H5C«S,HCCLS, helms, 
I H,H2,  H1I,MNN,N,HI,NP  1,!IP,SN,N!I0,MAXS,NS, 

I KHR,KHE,KHA. KAX,KBL,KEH,KHE,KXN, 

J BHS,KPL,KPU,KFF, 

> KNl,KMO,KRX,KCX,KXX,KFI, 

) KGB,KTX,KV1,KG1,KGA,KGB  ) 


Sabroutin«s  called: 


HOVE  MPS 


Called  by: 


HIMOS 


Purpose: 

To  input  constraint  data,  and  to  allocate  all  storage  except  that 

required  fcy  the  lu  file. 

Paraedters: 

BEA1*8  Z(HAXZ>  (Output)  The  constraint  aatrix  and  bounds  will 

he  in  varicus  par^s  of  Z,  pointed  to  by  the 
paraaeters  KEE ,KHA , KAX, KBL , KBU.  Certain  other 
array.”  will  also  contain  useful  inf oimation, 
caaelv,  those  pointed  to  by  KHB,KXN,K!iS. 

ISTEGBP  FEL,KPU,KPF,  ...,  KGA,KGB 

(Output)  These  are  starting  a'^dresses  in  Z for 
various  other  arrays  that  will  be  required  by 
EEIVEF  for  the  solvinq  procedures. 

Hethod : 

Three  calls  are  requited  tc  rubroutine  HPS.  For  the  first,  w;,  don't 

know  H or  V.  For  the  second  we  don't  know  K.  For  the  third,  the 

dlaensions  of  the  ptoblei  are  known. 

1.  Allocate  storage  it  Z for  rew  types  and  row  naees,  using  the  value 
HPCfcS  which  Bust  be  an  over-estlaate  of  the  actual  no.  of  rows  in 
the  data  about  tc  be  read. 

2.  Find  the  nearest  priae  nuabec  larger  than  MEOWS.  This  will  be  LBN, 
the  length  of  the  hash  table  fer  holding  row  naaes.  Allocate  that 
awount  of  storage  In  Z. 

4.  Allocate*  aore  storage  in  Z for  attays  hE,HA,A  (frua  the  bottca  up 
starting  ct  «rd  of  hash  table  KSTNAM)  , an.d  for  BL,  B0,1IAHF  1C,  S'AME?C 
(frcB  top  dowr.;.  At  this  stag*-,  El  and  BO  are  bounds  for  the 
slacks  only,  tc  be  defined  dutlrg  RANGES  input.  Also,  M'CLS  post 
be  ar  over-estiaat e of  the  no.  of  strurturals. 


Call  NFS  with  entry  2,  to  inpnt  COLUMNS,  BBS  and  RAHGBS. 

New  that  H and  N are  known,  wc  can  ccipact  storage.  First,  set  uu 
new  (and  final)  starting  addresses  for  existinq  arrays  and  for  sow* 
new  ones  HB,XN,H5  that  are  requited  during  BcnsDS  input.  ihen  copy 


existinq 


),XN,H5  that  are  requited  durirg  BcnsDS  input, 
krrays  into  their  new  positions  in  Z,  startin  i 


hen  copy 
fro  a the 


W 


bottOH.  In  particalnr,  the  esisting  bounds  on  slacks  nust  go  into 
the  end  of  tfce  nee  N-diBHr.sior.aI  arrays  PL  and  BU.  Word  alignient 
can  be  naintained  only  vith  the  help  of  subroutine  MOVE. 

7.  Call  BPS  with  entry  3,  to  input  BOUNDS  (on  structurals)  and  INITIAL 
EOUNDS  (on  nonlinear  variables). 

8.  Allocate  starting  addresses  in  Z for  all  other  arrays  required  by 
CFIVIP. 
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SOBBOOTINE  MBWPTC(EPS,  t , ETA,  ALAMTA,  U,  FU,  r,0, 

• XHIi,  FHIN,  f?!!!!*,  Xfc,  FH,  Gfc,  A,  E,  OLDF, 

• El,  SCXBE,  E,  t,  PP,  5S,  GTEST1,  GTEST2,  TOL, 

• ILOC,  IXPSl) 


Subrcntines  calltd:  Called  by: 

None  5EAPCH 


Purpose: 

To  find  a step  length  daring  a lineseatch.  For  full  docueentation  s*?'^ 
Gill,  Murray,  Ficken,  father  and  Wright,  NFL  Algorithms  Library,  HeE. 
No.  E4/16/F,  DNAC,  Natiocal  Physical  laboratory,  Teddingtor.  (Croirn 
Copyright  Reserved.) 


SCBPOUTINE  li15ECH(  S , I T I , IT 2,  »» MF 1 , K AH E2  , 

1 NCAPD,II0!:FND,?!A'<HSG,J1,J2,J?!A5|!,  J^CUKD  ) 

Subroutia«s  called:  Called  by: 

None  INSERT  LOADN  MPS 

Purpose: 

To  search  arrays  IPi  and  ID2  for  the  naae  defined  by  NAHE1  and  :vArB2. 
Naees  in  these  ▼ariafclet.  are  ir.  foraat. 


Paraaeters : 


INTEGER 

N 

(Input) 

INTEGER 

:di  fN) 

(Input) 

INTEGEF 

ID2  (N) 

(Input) 

INTEGER 

NAME1,N.ArJ2 

(Incut) 
locked  up 

INTEGER 

NCA8D 

(Input) 

INTEGER 

NOTFND 

(In ,out ) 
raae  was  : 

INTEGEF 

MAXHSG 

(Ir put) 
a Hewed . 

INTEGER 

J1 

(Input)  ! 

INTEGER 

.7  2 

(Input) 

T.N"’EGEF 

JMAFK 

(In, out)  1 
search  e 
search. 

INTEGER 

JPOOMD 

(Putput) 

Ctharwise 

No.  of  naaes  to  be  searched. 

Left  half  of  naaes  to  be  searched. 
Riqbt  half  of  naaes  to  be  searched. 
Left  and  right  halves  of  naae  hoinq 

• 

Card  no.  froa  which  NAHE1,NAHE2  came. 
Accuaulates  the  no.  of  tines  the  ir.pu*- 
net  foond. 

Lifit  on  no.  of  error  aet  sages 

HarKs  start  of  names  to  be  searched, 
Marks  end  of  names  to  he  searchod. 

On  input,  marks  where  the  pr^'vious 
nded.  On  output,  aarKs  end  of  current 

will  be  0 if  the  naae  was  not  found. 
, pcln’-.s  to  position  of  naae  in  IF  1,2. 


Method : 

Nothing  fancy  — just  a linear  search,  xplit  into  2 halves  (JHAFK  to 
.J2,  then  J1  to  JMASK)  , on  the  assuaptlon  that  the  naaes  being  looked 
up  are  likeiv  to  te  in  roughly  the  same  order  as  the  list  of  names  izi 
ID1  and  1!:2. 
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SOBPOOTINB  FACKIOI  H,M2,  »P,NAJ1 , •1AXU,IP1V, 

1 HP1?L,HPIVI,HIIVB,NI,NU,HL,HJJ,  OL,  Ct!,P,Y,IBGN,K  ) 


Subroutines  called: 


Called  by: 


None 


FACTOP 


Purpose: 

To  pack  up  alg  nificant  nonzeros  in  L,  0 and  F during  INVBRT.  The  ff-th 
ccluain  of  the  basis  has  just  been  processed  by  FACTOP  and  it  is  a 
spike. 


Paraneters: 


INTEGER 

IPIV 

(Input) 

Pivot  tow  for  the  spike. 

INTEGEB*2 

HPIVL  (N2) 

list  of 
will  be 

pivot  rows  for  the  columns 
added  to  the  list. 

of  L,  IFIV 

INTiGEa*2 

HPIVF  (H2) 

List  of  pivct  rows  for  F (previous 
lliv  will  be  added. 

spike  rows)  . 

INTFGERV2 

HPIVR  (*12) 

(Input) 

F ACTOH. 

The  full  row  permutation 

being  used  by 

RFALVfl 

Y (H) 

(Input) 

ronzetos 

The  transformed  spike,  containing  the 
tc  be  packed. 

INTEGER 

IBGN 

(Input) 

bump. 

Harks  ♦■he  beginning  of 

the  current 

INTEGER 

K 

Fcsition 

of  this  spike  in  the  basis. 

Method : 

1.  Increase  MSPIKE  by  1;  add  lElv  to  array  HFIVP.  Ose  HPIVF  to  pack 
up  nonzecos  in  all  spike  rows.  Reset  such  eleien«-s  T(i)  to  zero, 
so  they  will  not  be  packed  up  again  for  U.  Save  pivot  eleaent, 
?IV,  for  storing  in  the  L file. 

2.  Use  !1PI?R  to  pack  up  elesents  in  0.  These  are  in  roes  HPlVR(i), 
i=IPGV,  ...»  K. 

3.  There  cannot  be  any  eleaents  to  add  to  L if  this  is  the  last 

coluan  in  the  basis  (K=n) . Otherwise,  pack  remaining  elements  into 
1.  These  are  in  rows  HPIVP (i) , i=k*1,  ....  K.  Note:  The 

t rans*'ormation  is  going  to  be  used  as  if  the  nonzeros,  starting 
from  the  pivot  rcw,  were 

1.0,  Y(i1)/FIV,  Y(i2)/PIY, 

because  the  pivot  element  goes  into  P rather  than  L.  However, 
tc  avoid  the  divisicr.s  we  store  the  eleaents  directly  as 
PIf,  T(i1),  Y(i2), 

and  communicate  the  fact  tc  FTEANl,  BTPANL  by  setting  the  pivot 
index  HL(KL)  to  he  negative,  viz.  -IPIV. 
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S’JBEOr-TIllE  tPICEl  f,'!N,K,NF.,KF  1,S»i,NN0,Nr.,B»XS, 
1 JlA,HE,HS,HF,A,El,Bl),C,X,ri,GPC,G  ) 


Subroutines  called: 
Mod« 


called  by: 
cai VEF 


Purpose: 

To  select  ore  or  »ore  nonbasic  variables  eith  favorable  reduced 
oradieifs,  for  adoi*-lcr  tc  the  suuertasic  set.  The  class  of 
candidates  depends  on  cert  air  f ara^setars  as  felloes: 


1.  (PAKCAIO) 

2.  (fyLTIICH  PBIC^I 

3.  <r»FTIAL  Pi^ICE) 


If  ’<EEJ>C,  give  preference  tc  a list  of  variables 
in  array  JFE.1  that  were  rejected  frD®  the  basis 
cy  the  last  INVEFI. 

If.  K^ULPr’»C,  scar  all  non  basic  roluaas  and  select 
U].  to  S?1’JLFl<  candidates. 

Cthc;  wise,  scan  jest  ’•he  next  partition  of  A,  as 
well  as  ail  slacks. 


Faraveters: 

PiAL*3  PI{H) 

3PD{KH) 

IHT5f;EP*2 


?rAi*a  G(«iAXS) 


(Input)  The  pricing  vector. 

So*-  used. 

(In, cut)  If  K*L'lPP>0,  selected  cclunn  numbers 
will  be  add-ad  to  this  array,  starting  at 
HE  <'1S  + 1)  . 

(In, out)  If  ^r'JLFE>0,  the  cor  responding  r<  ducf-d 
grsdier. ts  will  le  added  to  this  array,  st, \rting 
at  G (NS+  !)  . 


>fet  hod : 

This  is  one  of  tho  mote  intricataly  (i.e.  , obscurely)  coded  routines, 
because  ’•he  various  options  ^ere  added  in  several  stages.  How'^ver, 


important  polr’ 

IS  that 

ir  all  cases,  the 

variables  J1,  JJ,  d3 

rsine  which  c 

ciuans  are  o be  priced  (s 

ea  the  loop  beiir.ning 

CO  dJ=j 1 ,G3  ) . 

11.'  cos 

sltie  values  are: 

J 1 

J2 

J3 

hawla.td 

1 

NSEJ 

NKEJ 

HULTIPLF  PPICE 

1 

NSIEDC 

J2  + M 

'-A5TIAL  PPTCE 

J1 

0 I’JiPHC 

J2tB 

where 


“RE.t  is  ’he  no.  of  rejected  variables  reaaini.ia; 

NST  EbC  = KP»?S”  1 warks  ’he  end  of  structurals; 

VFP'"  irTrgc/hrABFF  is  the  m.  of  columns  in  each  partition  cf  A. 

for  PIH'ilAI  PnICF  is  1t(a  multiple  of  NFFC) 
cr,  equivalently,  (previous  d2)*l  (mod  VSTHdO, 
for  the  last  two  cases  in  order  tc  include  the 
y s 1 :a  c k r . 


1.  If  PASTIAL  PPICB,  mcve  pointers  J1,J2  to  bracket  the  next  partition 
of  h. 

2.  Set  TOLD,  the  tolerance  to  be  ust'i  in  aeaturinq  whether  a reduced 

gradient  DJ  is  significant.  If  the  prohlfes  is  infeasible,  TCID  is 
always  TCLDJI  ♦PIHOF!'.  (Scaling  by  PlNCEri  is  net  rital  here,  except 
perhaps  if  a very  large  weight  is  used  *or  n®  ccuposit*^ 

objective.)  If  the  whole  ot  A is  going  to  be  scanned  during  this 

FBICE  (i.«.  if  VFiPPR^I  or  RULTIPIR  FSICE  is  in  effect),  we  sight 
as  well  allow  through  anything  that  is  going  to  be  regard  id  as 
significant  before  optiDality  is  declared.  Only  if  NPA?rK>l  is 
TOLD  going  to  be  decreased  d yaairically , by  setting  it  to  '"OLDJ 
which  is  reduced  syste  nintic  ally  as  desetiled  in  7 bHlow. 

3.  The  first  three  lines  in  loop  700  define  .1.  the  nex^  coluBt  to  be 
priced.  Skip  basics  and  superLasics  ittapd' ately . Skip  fixed 
nonbasics  nex*.  Then  treat  a coltiran  specially  if  it  is  i slack 
( DJ  = - appropriate  FI  value  ).  Otherwise,  th  - reduced  gradient  is 

DJ  * grad(.J)  - PT*A(d) 

where  grad(J)  is  0 or  C (J)  for  linear  and  ncilin.ear  variables  if 
the  solution  is  feasible,  cthotwise  0 for  .ill  v-,L-iabies, 

h.  Tha  required  sign  for  a DJ  depends  on  the  rttic  of  the  variable. 

Th'i  quantity  D is  defined  to  be  either  DJ  or  -DJ  or  .absfDJ),  such 
that  a positive  value  for  C moans  a f.avorable  t.’. 

5.  If  BOLTIFLB  PRICE,  eexpare  CJ  w^th  the  li.st  cf  vilues  in  G.  If 

necessary,  add  J and  DJ  to  HD  and  C,  .'-gch  ‘hat  ‘he  valuos  in  c 

re’Bain  in  descending  order  cf  sajni+’ude.  (This  section  was 
coded  by  Richard  Asmuth.) 

6.  Dpdata  CJSAX,  DJC,  JQ  tc  recerd  the  fees'-  DJ  found  so  far. 

7.  On  exit  fro»  loop  700,  detersine  if  any  suitable*  colunn  was  found 
(KONOPT>0).  If  so,  exit.  Otherwise,  "he  following  cases,  arise: 

3ANDAID  Eevert  tc  MGLTIFIE  or  oAhTXAL  3PICE, 

MULTIPLE  PRICE  Set  HEkSn  = numfeec  cf  columns  aided  ‘o  HP,  0; 

if  necessary,  delete  anv  small  entries  frcB  the 
end  of  these  lists. 

PARTIAL  PRICE  If  not  all  of  ^ has  beer,  seamed,  repeat  fro*  1. 

C‘herwise,  test  if  the  largest  DJ  ‘-•ncoun ‘t  red  in 
all  partitions  (CJSUP)  s smaller  than  the  level 
TC  inj  3*F  I>c  Eh , which  is  rakon  to  be  insior.  if  leant . 
If  sc,  exi*-.  ■Jtherwise,  r"'duct>  TGIDJ  to  a tenth 
ot  DJSUr  and  repeal  froir  l.  (At  least  one  colunn 

■ ust  get  picked  up  r.e/t  tvoe.  It  would  be  better 
1 to  exit  wirh  the  cclnmn  corresponding  o DJSUP 

I directlv.  Hewever,  the  reduction  of  TOLDJ  occurs 

only  orce  ct  twice  during  a tun.) 

I 

I 
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SUBROUTINE  PONCH(  H , N , KR HS , HS , BL , EU , Y , I D 1 , I r2  ) 


Sabcoutin«s  called:  Called  fcy: 

Nona  CBIVEF 

Purpose: 

To  ou^;  put  a basis  description  to  file  IfSCH,  in  a foraat  compatible 
with  soae  of  the  ccaaercial  LP  systems.  Each  carl  iaage  is  of  the 
fcra 


KEY  NAHE1 


nare; 


VALUE 


where  RBY  aay  be  IL,  CL,  XL,  XU  or  SB.  . .. 
of  the  standard  tcrsat,  necessary  for  description  of  norbasic 
solutions,) 


(Type  SB  is  a generalization 
for  description  < 

The  file  produced  by  EDNCH  nay  be  used  as  input  to  INSERT. 


Parameters: 
aEAl*3  Y (») 


INTEGER 

INTEGER 


ID1  (M) 
ID2  (N) 


(Input)  Contains  solution  values  for 
variables  (?t tuctotals,  RHS  and  logicais) . 
Vork  vector  to  hold  left  part  c£  row  names, 
hotk  vector  to  held  right  part  of  tow  names. 


all 


Method : 

1.  Read  row  names  fro*  the  scratch  file  Into  IDI,  1D2, 

2.  Process  each  structural  in  turn. 

(a)  Fead  colunn  name  from  scratch  file,  to  be  used  as  NABE1. 

(b)  (Nonbasics  or  supertasics)  Skip  nontasics  that  ate  at  a zero 
bound  or  are  fixed.  Cthetwise,  use  the  state  vector  HS  to 
index  the  local  array  NPY  to  get  the  correct  indicator  (LL,  UL 
or  SB).  Output  wi+h  NAME2  blank. 

(c)  (Basics)  Scan  the  state  vecter  tc  find  the  next  logical  that  is 
not  basic  (it  may  be  nontasic  or  superbasic) . The  corresponding 
tow  naae  will  be  NAHEi.  Select  the  correct  indicator  (XL  or  XU) 
«rd  output. 

Note:  XL  weans  the  logical  (not  the  row)  is  at  its  lower  bound. 

3.  Output  any  superbasic  Icgicals,  using  indicator  SB  and  NAR5.2  blank. 

It  IS  Inportant  that  these  fce  output  last. 

4.  Print  message  to  indicate  successful  PUNCH. 
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SOBROOTINR  P3  ( M ,n  1 , NB  E L 5 ?•  , KBU H T , HS  E IKE , HC,  HR , !) I , H J , 
1 HSC,HSR,HCLI£T,FBIIST,NSFIK'  ) 


E 

1 

i 

f 

t 

j. 


Subroutines  called:  Called  by: 

None  E4 

Purpose : 

TO  pereute  the  rows  and  ccluitr.s  of  a square  matrix  B so  that  it  is 
close  to  being  lower  triangular,  acart  trom  a fee  coluins  ("sril^es'M 
which  have  nonzeros  above  the  diagonal. 


Parameters; 


INTEGER 

E 

(Input)  Dimension  cf  B.  Will  be  the  size 
Must  one  of  the  bumps  in  the  whole  basis. 

of 

INTEGER 

PI1 

(Incut)  f1*1. 

INTEGEF 

NBELEM 

(Input)  No.  cf  nonzatos  in  P. 

INTEGER 

N3iJf!P 

(Input)  The  r.uitbet  of  this  bump,  used  only 
storing  into  array  HSPIKE. 

for 

INTEGEP*2 

PSPIK'  (M) 

(Output)  HSfrKt(t)  will  he  -NPiJNP  if  column 
is  a spike,  ♦NuUMP  ctherwise. 

1 

INTEGER*? 

HC  (B) 

(Output)  Tie  assigned  column  order. 

INTEGER*? 

HB  (H) 

(Output)  The  assigned  row  order. 

INTEGER*? 

HI(B) 

(Input)  Ccluin  counts. 

INTEGER*? 

HJ  (H) 

Work  vector  to  hold  row  counts. 

INTEGER*? 

H5C  (M  1) 

(Input)  HSC(-j)  points  to  the  start  of  column  j 
in  HCLI3T. 

INTEGER*? 

HSR  (H  1) 

work  vfitor  to  point  to  start  of  rows 
(tRLIST. 

in 

INTEGER*? 

HCLIST  (NBEL' 

“) 

(Inuut)  Cclunn  list. 

INTEGEF*? 

HRLIST (NPEIE 

Hcrk  vector  to  held  tow  list. 
INTEGEB  KSPIK’='  iCutput)  Fo.  of  spikes  found. 


Method : 


The  bulk  of  P3  was  coded  directly  frem  the  description  of  the 
Preassigned  Pivot  Frccedurc  by  Hellerman  and  Barick,  Math. 
Prograeming  1,  2,  November  1971.  It  has  t*=en  eodified  as  follows: 

1.  The  backward  and  forward  triangle  sections  have  been  deleted,  since 
they  are  now  handled  by  P4.  nance,  this  version  of  P3  is  intended 
for  use  on  an  irreducible  square  watrix. 
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SOBBOOTI«E  i?4  ( K , H 1 , H,  ME , HP  1 , B Z2  , HSPIKE  , !IA  , HF  , HS  HP  , H I , J, 

1 llUHB,HP,B,1i1,M2,HFM,NEUHE,MSPIKE  ) 


Subroutines  called:  Called  by: 

BOHPS  HKLIST  P3  IPNSVl  IMVEPT 


Purpose: 

< To  paZBuce  the  roas  and  coluant-  of  the  basis  eatrlx  B into  the  bunp 

! and  spike  crdaring  of  Helletean  and  Barlck.  The  perauted  B Mill  be 

block  lower  triangular. 


Paraaeters: 

IHTEGEfi  MZ2 
INTEGEF*2  HSPIKE(H) 


COO  1 1-1-10022-2-2000 

which  Hculd  aecn  that  the  last  two  coluans  ia 
each  bump  are  spikes,  the  forward  and  backward 
triangles  are  both  of  size  1,  and  the  two  bueos 
are  separated  by  2 triangle  coluans. 

HITEGEP*2  HB  (H)  Cc  input,  contair.s  column  nos.  of  basic  wars. 

INVEFT  sets  this  list  up  backwards,  so  slacks 
appear  first. 

IWTEGEFm2  BC  (H*  (Cutput)  Will  contain  the  coluen  nos.  in 

pereuted  order,  ready  for  FACTOR. 

HiTEGEP*2  HR  (fl)  (Cutput)  will  specify  the  preassigr.ed  pivot 

tows  for  the  tacrdered  columns. 

IBTEGER*2  HI  (H. ) , H J (H ) , NUHB  (H ) , HE  (M ) , E (H ) , H 1 (M)  , W 2 (H) 

All  wcrk  vectors. 

I1ITEGE3*2  HFN(NZ2)  Work  vector,  will  be  used  for  coluan  and  tow 

lists  by  P3. 

IMTEGEP  NBOHP  (Output)  Final  nc.  of  buaps. 

IITEGSP  NSPIKE  (Cutput)  Final  no.  of  spikes  in  all  buaps. 


(Input)  Twice  the  no.  of  nonzeros  in  E. 
(Cutput)  Codes  the  busp  and  spike  structure  in 
the  following  way.  Suppose  B turns  out  to 
have  2 buwps  of  dirension  4,  each  with  2 
spikes.  HSPI KE  sight  then  be 


'i 

ij 


Ret hod: 

1.  Set  up  coluen  list  in  order  specified  by  HP.  This  order  is 
arbitrary  exceot  for  safety  the  slacks  ace  listed  first. 

2.  Call  TBMSVL  to  find  a eaxieal  transvecscl.  HP  will  cortair.  the 
resulting  column  pertutation. 

3.  If  P ia  structurally  singular  (HHEJ. ME. 0) , an  indirect  lethod  is 
used  to  find  which  columns  in  HE  were  selected  and  which  of  those 
Bust  be  replaced  by  suitable  slacks. 

4.  Reorder  HB  and  make  a new  coluan  list. 
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5.  Call  BUBFS  to  find  a symetiic  p^^rautation  HP  that  aakps  B block 
upper  triangular. 

6.  Since  we  want  B to  he  block  lower  triangular,  we  now  process  the 
output  froa  B0(*PS  in  teverfc  order,  in  particular  the  array 
B(NBLK)  which  points  to  the  begianirq  of  each  huap  within  array 
HE. 

7.  For  each  buap,  pointers  Fl,F2  aark  the  beginning  and  end  going 
backwards,  and  Q0,C1  aark  the  beginning  and  end  going  forwards  when 
results  are  stored  finally  in  HC  and  HB. 

8.  Buaps  of  diaension  1 or  2 are  treated  directly.  Otherwise,  the 
ccluan  list  for  a bump  is  set  up  by  taking  the  tow  nos.  in  each 
column  in  the  bump  and  csing  thea  to  index  the  array  HUHD  which  is 
output  by  BUMPS.  If  the  relevant  value  of  NUMB  lies  between  Pi  and 
P2  then  the  element  lies  inside  the  butp. 

9.  Call  P3  and  accuaulate  results  forwards  In  HSPikE,  HC,  HB. 


SUEBOOTIME  BESETF  ( N M , NP  ,IIS,  P , COND  ) 

Subccutia«s  call«d:  Called  by: 

None  rSIVEP  invert  BRITN 

Purpose: 

Besets  E to  the  identity  latrii  (B  upper  triangular*  stored  by 
coluins).  Called  when  the  first  feasible  point  is  found,  or  as  one  of 
the  recotery  actions  if  the  lineseerch  fails. 

Paraieters: 

BIAL*P  COMP  (Output)  Beset  tc  1.0. 
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SOBFOOTIIIE  FGITSf  P! , «2  , M S , (INN  , N , NE,  NP 1 , HP,  Nll»  ■»  0,  H*XS  , NS  , 

1 -AXi,,BAXr), 

2 HSPlKE,HA,HB,FE,HS,HFIVL,HPIV0,HPI7F,NL,NU,A,il,B0,i,PI,X, 
J GEC,Y,Y  1,XII,C,G1,G,CNFW, 

4 Hl,ni.,Hn,RD,  F ) 


Subccutines  called: 


Called  by: 


ADDCOL  AlIGN  BTFANU  CALCG 
CG  CHGZg  CHtrZR  CCPDPF 

ORLCCL  ECTI  PTPANI  FUNGBD 
GETGFD  BiJSSIR  BTRSOL  B 1V30D 
SEARCH  SETPI  SQHEEZ  UNPACK 


Purpcse : 


tSIVEF 


To  pcrfot*  an  iteration  of  the  reduced-gradient  aethod. 


Parareters: 


All  defined  previously. 


Kethcd : 


SPHS  is  either  3 or  4 on  entry  to  EGITN.  The  only  difference  is  that 
when  NPHS=3,  a list  of  nontasic  variables  (set  up  by  PBICE)  aast  first 
be  added  to  the  supertaric  set. 

Thereafter  (for  both  phases)  tie  routine  perforas  one  optioiization 
step  on  the  curr-at  basics  and  superlrasics,  la  general  i'  aay  be 
entered  with  NPHS=4  for  several  consecutive  iterations,  during  which 
the  nuvber  of  sucarhasics  will  either  regain  constant  or  decrees--. 

1,  (Initialize)  Set  up  various  Icgical  variables  to  deterviie  the 
current  state  of  the  probleir  (fEASBl,  INFSP.L,  LINEAR,  NDNLTN, 
VAPnET).  For  exaiple,  LIN'AP  is  true  if  the  current  solution  is 
infeasible,  or  i*  there  are  no  nonlinear  variaMas;  the  linesearoh 
and  guarl-Newton  or  conjugate-gradient  steps  can  then  be  stippe.d. 


(Add  superbasics)  If  NPH 
H3(j)  , 

beccae  superfcasic.  The 
reduced  gradient,  TJC. 
BGKOPf!  (the  nor*  cf  th 
super  basics)  , Ignore  t 
Otherwise,  add  the  ne 
accordir.  dlj.  If  neces.<3a 
exit  If  NS  has  alrea 
exceeds  NAIF  (ftessian  di 


£=3,  the  nontasic  variablas  listed  in 
KtNStSEWSP  have  been  set  up  by  PPICE  to 
first  variable  in  the  list  has  the  largest 
if  DJQ  is  not  sufficiently  l«rq«t  than 
o reduced-oradlent  vector  for  the  existing 
he  list  and  proceed  as  if  NPHS=4, 
w superlasics  one  by  one  and  incr*  amt  NS 
ry , add  new  coluans  to  R.  Take  failure 
dy  reached  NAXS.  Set  VARHTT  = EAIS-:  if  NS 
oersion) . 


(Coapute  search  directicn) 

(a)  Generate  a search  direction  s for  the  snperbaslcs,  storinq  s 
in  Y (H«-i)  , j-  1,  ....  (iS. 

If  lINEAB,  s Is  the  negative  reduced  gradient  vector,  -G. 

Tf  fAlEET,  f<*fs  ' -G. 


no 


If  COHJGB,  s will  already  have  teen  set  up  by  the  previ  as 
entry  to  KGITS,  except  if  PISTPT  = TfUE. 


h- 


i 


It 

It 


(b)  Solve  By  * -Sfi  tc  yat  cc  r res  pcndinq  secirch  direction  y for 
basics.  Store  y in  Id)*  j-1»  ...»  M. 

<J.  (CHU7.S)  Find  which  variable  teaches  one  cf  its  bounds  first  when  JC 

is  changed  to  X < TtS7A*Y  (TtFTA  >=  0).  Set  JP  accordingly, 

treble*  is  urtounded  it  iTNEAt  and  .7P=0. 

5.  (Linesearch)  Sfcip  if  LJNEAE  or  if  THETA  is  essentially  zoro. 

Otherwise*  call  SiABCH  tc  detetsine  a new  step  length  THETA. 

Declare  problea  untoiinded  if  recassacy.  If  SFAPCf’  fails  on  one  or 
sore  consecutive  iteratiers,  try  the  following  in  turn: 

(a)  Reset  the  Hessian  (B  - I)  and  try  again. 

(b)  Call  INVEST  and  try  again. 

(c)  Switch  to  VPH£=3  tc  ask  for  PRICE  (sore  superbasics)  and  try 
again. 

(d)  Repeat  (c)  up  to  u tlses. 

(a)  Give  up  for  ‘■he  day.  Try  again  teserrow. 

6.  (Find  new  reduced  gradient)  Skip  if  LINEAR  or  THETA=0.  Otherwise, 
cowpute  new  PI  vecter  and  new  reduced-gradient  vector  GNEW.  If 
VAFfiET,  coopene  T = GNEH  - C,  etc.  rea'ty  for  quasi-  Newten  update 
to  Hessian. 

7.  (Basis  change)  If  i <=  jp  <=  a basic  variable  reached  its 

bound. 

(a)  If  VABH37,  do  quar.i-Newtcr.  update  tc  B. 

(D)  Solve  P’y  = e (JF)  . 

(c)  Call  CHUTC  tc  find  the  Kyth  supertasic  to  replace  basic. 

This  will  he  rclv^n  no.  JvK 

(d)  Tf  call  “ P5CD  tc  aodify  f to  account  for  basis  change. 

(e)  dodify  PI  iisirg  y fees  (b)  above. 

(f)  Solve  L*Y1  - A(Jv'1)  to  prepare  for  update  of  LU  factors 
of  basis. 

(g)  Go  to  3 (b)  . 

8.  (Delete  one  superbasic)  It  JP  > B*  a superbasic  variable  reached 
its  bound.  Set  j<Q  * gg  - g, 

(a)  If  VAKHET,  do  quar t-Newton  update  to  E. 

<b)  Cali  DELCri,  to  delete  the  KQ-th  superbasic  froa  X,  H"  and 
perhaps  K. 

(c)  Becoapute  re iuced-gradaent  G ard  nora  FGhOPM. 

)d)  Go  to  17. 

9.  ('jneonstraine d step)  If  JF  = C,  a •ininua  was  found  along  X ♦ 

THfTA*T  before  anv  variable  hit  a hourd. 

(a)  If  COMJGE,  cill  CG  to  confute  search  direction  for  superbasics 

tor  next  | 

<b)  If  VARHEt,  do  quasi-hewter  update  to  F.  il 

(c)  Set  G •=^  GNFii.  ^ 

10.  fEstiaate  condition  no.  cf  reduced  Hessian)  Ccapute  DHIX  and  DHIN* 

the  larqe.st  and  srailest  dlagcral  eleaer.ts  cf  P.  Set  COND  = 

(Dr AX/DBIN) •*?. 

I 

11.  (Convergence  test  for  current  solspare)  Deterfeine  whether  to  stay 

in  phase  **  (ano^  'er  iteration  with  the  saae  suoerbhsics)  ct  to  as*  • 
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SOBFODTIHi:  HTRSCK  P,y,NF,N  ) 


Subroutlo«8  called: 

None 

Purpose: 

To  sol Ye  the  systea  B'Ey  » z, 


Called  by: 
PGITN 


Paraietcrs: 

INT3GBH  I! 

K2AL*8  r(V) 


Method : 


(Input)  Current  dimension  of  B. 
Initially  holds  2.  overwritten  by 


Forward  substitution  followed  ty  bacK-substitutlon.  Both  are  arranged 
to  run  through  the  cclumas  of  B,  in  order  to  access  meaorv 
sequentially.  This  iaptcves  execution  speed  on  sone  aachines,  e.g. 
those  with  a cache  memcry. 


SOBFOOTINE  PUDC(  F , ¥ , NF ,11 , TOLZ  ) 


Sabroatines  called:  Called  by: 

None  FI  ROD  FI  PROD 


Purpose: 

To  aodify  R such  that  (nen  B*B)  = (old  P'B)  ♦ Yy*» 


Paraaaters: 

INTEGSB  R (Input)  Current  dimension  of  R. 

REAL*8  T(N)  Cn  input,  contains  vector  y.  Gets  overwritten. 

REAL*8  ICLT  (Input)  Tolerance  for  skitpinq  transfer nations. 


Method : 

A single  forward  sweep  of  plane  rotations.  For  exanple,  ¥(1)  is 
eliiinated  first,  changing  the  renainder  of  T and  th#>  first  row  of  P. 
Maintains  nonsingularity  cf  R,  Undt»rflcw  way  occur  occassionallv 
witheut  fatal  consequence. 
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SUDBOaTINE  P1MOD(  F , Y, NB ,N , lOl 2 , SIGMA , K ) 


Subroatines  called: 
BUCD  P1SUB 


Calltid  by: 


COHrFF 


Purpose: 

To  Kodify  B such  that  (new  R'B)  - (old  R*P)  ♦ SI6!!A*yy*. 


P^raiceters: 


INTEGER 

RFAL*3 

RFAL*3 

FEAL*8 

INTEGER 


M 

T (N) 
TOLZ 
SIGMA 

K 


(Icfut)  Current  dimension  of  P. 

Cn  input#  contains  vector  y.  Gets  overwritten. 
(Input)  Tolerance  for  skipping  transformations. 
(Input)  A scalar,  may  be  positive  or  negative. 
(Output)  Ettor  flag,  set  to  -K  if  modification 
f aile  d. 


Method : 

1.  Change  y to  abs  (SIGMA)  ty.  Exit  if  non  (y)  <=  TOLZ. 

2.  Call  PIADD  or  P1SU8  according  to  the  sign  of  SIGMA. 
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SOPPODTIIE  aip«cr(  F,V,Y,NB,H,'I01Z,T0LD,  JQ.KADD  ) 


Subicutines  call€d:  Called  by: 

P1KDC  FGITN 


Purpose : 

To  Bodify  H such  that  (n«?«  B'f)  » (I  ♦ TB*)*(old  ♦ wv'J,  for 

the  special  case  where  w is  qeneral  I ut  w is  the  JC-th  unit 
yectcr,  e (JQ)  . 


Paraietecs: 


INTEGSB 

N 

(Input)  Current  diaension  of  R. 

• r 

REAL*8 

Y(») 

Cc  input,  con*'ains  vector  y.  Gets  overwritten. 

PTAL*ft 

TOLZ 

(Input)  Tolerance  for  skipping  transfer ma* ions. 

BEAL*d 

TOLD 

(Input)  Miniaua  allowable  size  of  the  J0“th 
diagonal  of  f. 

ISTEGEB 

FADD 

(Cutput)  Error  flag,  set  to  -KADD  if  the  JQ-th 
diagonal  of  F had  to  be  raised  to  TOLD. 

Method : 

1.  JC«1  is  a special  case.  siaply  add  a(1,1)*v'  to  the  first  row  of 
B and  exit. 

2.  Otherwise,  compute  y = Fw  * Be (JQ) . i.e.  pull  out  the  JQ-th 

ccluan  of  B. 

3.  Perfoi*  a partial  backward  sweep  of  plane  rotations  (as  in  FISOB) 

reducing  the  vector  ( y 0 ) to  { 0 s ) by  rotations  in 

the  clanes  ...»  1.  This  creates  a row  r 

below  P. 

h.  Coapute  y * r ♦ Sv. 

5.  Give  y to  BIAPC  to  coaplete  the  aodi f Ication. 

6.  Test  the  JQ-th  diaucnal  of  F,  which  could  he  rero  or  ?ery  srall. 


Postscript:  This  routine  has  since  been  modified  to  perform  a partial  backward 
sweep  and  then  just  a partial  forward  sweep  (instead  of  the  caii 
to  RIADD) . See  Murtash  and  Saunders  (1978)  . 


F 
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SUBHOUTTN?,  B1SUE(  B , T , NF  ,»  , TOIZ  , F *I  L ) 


Subroutines  called:  Called  by: 

None  FI  HOD 


Purpose; 

To  Bodify  P such  that  |ne»i  P'F)  = (old  R*P)  - yy'. 


Peranaters: 


INTKGSE  n 

SFAI*8  Y(N) 

REAL*8  TOLZ 

LOGICAL  FAIL 


(Input)  Current  diacnslon  of  B. 

Cn  input,  contains  vector  y.  Gets  overwritten. 
(Input)  Tcierance  for  skipping  transfonoatiors. 
(Cutput)  Pill  be  TP13E  if  y is  such  that  (new 
F’P)  would  not  be  positive  definite;  in  this 
case,  P will  te  unaltec*^"!. 


Method ; 

1.  Solve  B'p  = y.  Coapute  D * 1 - p'p, 

2.  Ejit  with  FAIL  = TPCF  if  D <=  TCIZ.  Otherwise  set  S = SORT  (D) . 

3.  Ferfora  backward  sweep  of  plane  rotations,  reducing  the  vector 

f F S ) to  ( 0 1 ) by  rotations  in  the  planes  (i,N*1), 

:=N,N-1,  ...,  1. 


j 

1 

1 
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SaEPOOTIliB  SAVES  ( K , H , NNO  , MS,  «S , HS,  HE,  X , XN  , ISTATE  ) 


Subroutines  called:  Called  by: 

Nona  CFIVEP 


Purpose: 

To  output  the  current  state  of  each  variable  in  the  form  of  a bit  map 
(one  character  per  variable)  , along  vith  coluan  nutibers  and  values  for 
supetbasic  variables.  Output  is  to  file  IBEKB. 


ParaRdters: 

IHTEGEP  ISTATE(3)  (Input)  A 12-character  string  dtscribinq  the 

current  state  of  the  solution  fsee  subroutine 
STATE) . 


Method : 

1.  Output  one  card  image  containing  FrobleR  naae.  Iteration  number, 
ISTATE,  Phase,  Oh;^ectiv€  value  (cr  Sum  of  iafeaslbilities,  if 
Phase=  1)  . 

2.  Output  one  card  inaoe  containino  OBJ  naae,  ?HS  naae,  BANGv  rane, 
PC'JND  naae.  No.  cf  ccnstraints  (P)  , No.  cf  variaLiss  (N)  , ho.  of 
superbasics  (NS) - 

3.  Output  the  state  vector  US  using  FOPHAl  (801  1)  . 

4.  Output  the  column  nuater  HP  (J)  and  value  X (J)  for  each  superbasic, 
cne  pair  oer  card  i-rage,  using  FC’=MRT(1P,  1P224.  14).  Terainato  the 
list  aith  coluan  number  0. 

5.  Peaind  file  INEWE. 

f».  Print  message  ‘BASIS  PAE  SAVED  ON  FILE  <INE«B>' . 
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SUPFCDTHI2  SEJFCHf  NS  , !1S,F  N , X , P ,GBD  , HB  , X K ,C  , X 1 ,(3  1 , 
1 ETA^PNOBf.STI EHX,ALEH«.F,UNBDC,IFAIL  ) 


Subroutines  callad:  Called  by: 

DOT  FTNGFD  GETGPC  SEWETC  FGITS 


I Purpose: 

■ I 

; I To  find  an  approximaticn  tc  the  miniaui  (sith  respect  to  ALPHA)  of  the 

; function  f (X  + ALPBA*P)  defined  by  FUNGSD.  Driver  for  NtwpTC, 

derived  directly  trci  subroutine  LNSRCH  by  Gill,  Hurray,  FicKer, 
Barter  and  Wright,  KFl  Algctithas  Library,  Pef.  No.  E«</16/r,  DNAC, 
Nationa]  Physical  Laboratory,  Teddington. 


Paranetars 

PFAl*fl 

X(HS) 

BFAL*8 

P (MS) 

“FAl*^ 

X 1 (HS) 

BFAl*8 

G 1 (NN) 

REALMS 

ETA 

HF Ai^a 

PNOfiH 

REAL*8 

SIEPHX 

RFAL*8 

ALPHA 

PEAL*8 

F 

LCCTC AL 

ONBDD 

INTEGER 

IFAIL 

On  input,  contains  initial  values  of  the  basic 
and  superbasic  variables.  On  output,  contains 
X * ALFHA+P. 

(Input)  Contains  the  search  direction  for 
basics  and  suparbasics. 

Work  vector  used  to  hold  teanorary  X ALi  iA*P. 
Work  vector  to  held  gradient  for  nonlinears  at 
X 1 . 

(Input)  Tbe  linesearch  accuracy  tolerance. 
(Input)  Nora  of  P. 

(Input)  Upper  liaiit  on  ALPHA,  defined  by  the 
bounds  cn  X. 

(Output)  Beturns  step  length. 

(Cutput)  Feturns  best  ob-jective  value, 

(Output)  Will  be  TFDE  if  ALPHA  turns  out  to  be 
very  large  ( I.OFefi  currently  used  ). 

(Cutput)  will  ne  0 if  search  was  successful, 
positive  otherwise. 
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SCBBOUTINB  SETPl  ( f , H k , N E , KE  1 , 8P , « A XL, ftA XU , NORM , TOLZ, 

1 HSPIKE,HA,HE ,FP1VI , HE  I VD , HF IVP, N I , N U, A, P I, Y , 

2 HL,DL,H'J,DU,  F ) 


Subroutines  called: 


BTFANL  PTFAWa 


Called  by: 


ADDCOl  DRIVER  RGITH 


Purpose : 

To  solve  the  systes  E (transrcse) •?!  * Y,  given  the  vector  Y, 


Paraieters: 

IBTEGEF  NOBS 

aEAI*B  TOlZ 


PSAL*8 


HEAl*8 


Hethcd: 


FI(H) 


(Input)  0 if  norm  (PI)  is  not  required. 

(Input)  If  r.otB(PI)  is  coaputei,  any  PI(i) 
smaller  than  TOLZ  mill  fca  reset  to  zero. 

(Output)  The  required  solution.  (The  actual 
parameter  may  net  always  be  the  pricing  vector 
for  the  sinplex  method.) 

(Input)  The  given  r)i£-  vector.  No<-  overwritten. 


actual 


1.  Set  PI  = T and  call  PTPAhO,  BTRAKI. 

2.  If  SORS  is  nonzero,  compute  SDK  = sub  of  PI(i)**2,  resetting  saall 
components  to  zero  and  counting  significant  components.  Set  riMOR" 
= square  root  cf  the  avecaoe  value  of  nonzero  ?I(i)**2.  This  is  an 
attempt  dampen  the  effect,  of  a few  very  large  components. 


SOBFOUTIHE  oEIX(  H ,!I2  , M,  IF  , MF  1 , IF  , NS  , HS  , MiXL,  (I  AZU, 

1 HA,BB,HP, HS, HFIVl,Hf 1?0,HP1»F,NL,»D,A ,PL,BO,X,T,I1, 

2 HL,l)I,4U,nU,F,IFESET  ) 


Subroutines  called: 


BTPANU  FTFASl 


Called  by: 


CRIVEF  INVEPT 


Purpose : 

(Optionally)  to  ccapute  <ialu€s  for  the  basic  variables^  using  the 
current  fcasia  factorization  and  values  for  superbasics.  Also,  to 
perfctn  a row  check  cr  X,  to  ensure  that  the  general  constraints  Ax 
= b ire  satisfied  to  a certain  precisicn. 


Faraaeters: 


EEAL^B  X(MS) 


EEAL*«  T (II) 

BEAI^P  T1(*5) 

INTEGER  T.RESFT 


Supertasic  values  are  input  in  X(1),  1=Hel, 
?*?.  Values  for  basic  variables  are  either 
input  in  X ( j)  , j = ''»  •••»  M»  or  are  to  be 

cctputed,  depending  on  the  value  of  IPE3ET. 
lork  vector. 

Work  vector. 

(Input)  0 it  basic  X(j)  are  input,  >0  if  they 
ate  to  be  reccaputed. 


Kethod : 

1.  Set  Y = 0. 

2.  For  each  nonbasic  coluvn  A (F) , accuiulate  -A{K)*B  into  Y, 

where  B is  the  acptcpciate  upper  or  lower  bound,  depending  on  the 
state  ot  the  associated  variable.  Can  he  skipped  if  B is  zero  (as 
it  usually  will  be) . 

3.  Fcr  each  supartasic  cclunn  A (K)  , accuaulate  -A(K)*X(J)  into  T, 
where  X (J)  is  the  aporopriate  value  of  the  supertasic  variable, 

4.  (Fecoapute  X)  If  TF?SE1>0,  solve  B*X  = Y for  the  basic  vars. 

5.  (Iterative  refineaent)  If  IBE3ET>0,  set  ¥i  = Y - B*X  and  solve 

the  systea  P’CX  * Y1  fcr  a correction  to  X.  Store  PX  in  Y1. 
Ccapute  DXHOF'"  = the  largest  coapontnt  of  DX,  excluding 

corrections  to  Icgicals. 

6.  (Fow  check)  Set  Y * Y - P*l,  At  the  saae  tiae,  ccapute  ZNOPW  = 

the  largest  coaponent  of  \,  excluding  loqicale;  also  compute  Yl(i) 

= ncra  of  the  i*tb  row  of  E. 

7.  If  iterative  reflnetent  was  done,  an  astiaate  of  the  condition 
ruaher  of  B is  r XICFF/ (XNCBH^SPS) . 

B,  Ccapute  the  largest  normalized  rcw  residual,  BHAX  = 

T(i)/Yl(i).  Set  the  error  flag  IEPP*6  if  FflAX/XNOSn  is  larger 
than  the  row-check  tcletarce,  ICLIOW. 


SUBBOUTINB  SCLH  ( SCAll  ,M  ,N , NB,  NE  1, NH ,!»|10,  Mb,  KS, 

1 A,E1,BU,XN,C,X,EI,T,HST,1STATE,0IIMSK  ) 


Subroutines  called; 


SCLPBT  DMPACK 


Called  bv: 


CRIVEF 


Purpose ; 

To  load  solution  and  state  values  for  all  variables  into  arrays  Y and 
HST,  and  to  output  the  soluticr,  in  apptcxiaately  BPS  foraat. 


Parameters: 


IKTEGEP  NCAIL 


PIA1*B  Y{N) 


IMTEGSB*2  RST(M) 


(Input)  NCALl=1  if  Y and  RST  are  to  be  set  up; 

NCALL=2  if  solution  is  to  be  output. 
Rclds  solution  values  as  follows: 

Y(1),  ....  Y(KPHS-I)  Strocturals 

y(KRHS)  FHS,  always  -1.0 

Y(KFHStl),  1 (M)  Logicals  (not  rows) 

Rclds  state  values  for  structurals,  PHS  ai.? 
Icgicals,  as  fellows: 


State 


IRTEGEP  IST»TE(3) 

LOGICAL  CMDISK 


Weaning 

Variable  is  at  its  lower  limit 

Variable  is  at  its  upper  liait 

Supertasic 

Basic 

Fixed 

Free  (but  nonfcasic) 

Below  lower  bound 
Above  upper  bound 


(Input)  Contains  a string  describing  the  state 
cf  the  solution. 

(Input)  TFOE  if  solution  is  to  be  output  ^o 
file  ISCLN. 


Rethod : 

1.  (MCILL^I) 

(a)  Copy  the  state  vector  RS  into  array  HST. 

(t)  Pun  through  the  ncnhasics,  setting  Y to  the  appropriate  bound 
values  and  changing  LL  or  UL  states  to  EQ  or  FR. 

(c)  Pun  through  the  basics  and  supertasics,  setting  T to  the 
appropriate  X values  and  changing  ES  or  SBS  states  to  --  or  ♦♦ 
for  any  infeasible  values. 

(d)  Exit  if  ptoblea  is  linear  or  infeasiMe*  or  if  usee  has  not 
requested  one  last  call  bo  CiLCFG.  Otherwise,  overwrite  any 
nonlinear  scluticn  values  with  the  values  in  XM,  in  case  the 
user  altered  XM. 
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2.  (NClkLL>1) 

(a)  Set  LPE  to  thf  the  output  unit  (file  t or  filfa  ISCLN). 

(h)  Expand  the  RHS  vector  into  X. 

(c)  Output  solution  heading  (NAME,  OBJ  value,  STATE,  PHASE,  etc.). 

(d)  Output  ROMS  heading. 

(e)  For  each  row,  convert  state  and  solution  values  for  logicals 
into  those  for  rows  (Rows  valua  * FHS  value  - Logical  value); 
read  the  Row  Nane  froa  file  ISCB;  rail  SOLPPT  to  output  me 
line  of  the  solution. 

(f)  Output  COLDHMS  heading, 

(g)  For  each  colusn  J,  ccs.pute  the  reduced  gradient  DJ  by  the  same 
method  as  used  in  £?1CE;  include  the  gradient  value  C (J)  if 
solution  is  faasitle  and  variable  J is  ncnlinear;  change  any 
negligible  values  to  zero;  read  the  Column  Naae  from  file 
ISCR;  call  30LFFT  to  cutput  cne  line  of  the  solution. 
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SUBBCOTISE  SCIPS1(  .1 , 1 C , lA  , X , C . B 1 , B 2 , D , K , P f lUS, C NDISK  ) 


Sttbioutines  called:  Called  by: 

Non#  SOIM 


Put  pcs« : 

To  output  ona  line  cf  the  sclution.  Used  for  both  the  BOWS  and  the 
CCLUB.NS  sections. 


Paraaeters:  (All  are  input  pacaeeters) 


IK’"EC-!iP  J 


INTEGER  ID  (2) 

integer  ia 

REAL*8  X 

REAL*8  C 

BEAL  E1 

HEAL  P2 

FFAL*(J  C 

INTEGBB  K 


PEAL  BPL'JS 

LOGICAL  CNDISK 


The  variable  nuabet. 

for  the  first  structural; 

J=N  tor  the  last  tow. 

NaBie  cf  tb*-'  tew  or  ccluan. 

An  incex  defitinq  the  state  of  the  variable, 

1 he  value  of  the  variable. 

Elar^  activity  (for  rows)  or  objective  gradient 
(for  colutr. s)  . 

Lower  found  on  row  or  colnan. 

Ur  pc  r tour.d  or.  r rw  or  colu  an . 

Fi  value  (for  rows)  or  reduced  gradi-^at  (for 
columns)  . 

Ihe  ccluan  no.  used  by  HPS/360. 

K=1  fer  the  first  row; 

K=N  ter  the  last  structural  (the  FHS) . 

"Plus  infinity'*  for  I curds. 

TFUE  if  the  scluticn  is  being  output  to  file 
ISCLM,  rather  than  the*  printer. 


Be*-  hod; 

The  only  logic  here  is  in  detenining  ahich  format  stateaent  to  us>, 
depending  on  the  size  of  the  two  bound  values,  El  and  B2.  If  CNDISK 
is  IKUF,  output  occurs  directly  to  file  ISCLN,  using  IPElfc.h  for  ail 
BEAL  quantities.  Ctharwisn,  cne  of  focr  possible  format  stat  ^aents  is 
used.  The  word  "NCNE"  is  used  for  all  infinite  bounds;  otherwise, 
BEALG  are  output  with  foraat  P16.h. 
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SOEBOOTINE  SPECS  ( I , I SPEC S , ISCF , I VEBSN , NSPRC  ) 


Subroutines  called:  Called  by: 

Hon#  flINO? 


Purpose: 

To  Bake  a first  pass  at  processing  the  SPECS  filt-.  Each  card  of  ♦he 
file  is  assuied  to  contain  one  of  the  fcllowing: 


KEIWOBD 

EZYBOPD 

KEYWORD 


IDERTIFISB 

IMIEGEE 

REAL 


in  rtJlatiTcly  free  foraat.  ihe  aiia  here  is  to  extract  such  iteas  as 
certain  character  strings  and  to  output  them  to  the  scratch  tile  ir  a 
fixed  fcriat  tha*-  is  suitable  fet  re-reat’ing  using  standard  I/O 
routines.  This  is  cna  vav  of  processing  da*^a  in  free  foriat,  in  a way 
that  is  Machine-independent. 


It  is  assuaed  that  the  first  and  last  cards  in  the  SPECS  file  contain 
the  KEYWOFDS  "BERIh'*  and  "ENE"  respectively. 


I Paraaetecs: 

I INTEGER  1(80) 

f INTEGEB  ISPECS 

I INTEGER  ISCF 

\ INTEGER  I?EBSN 

t INTEGER  NSPEC 


Korkspace  to  hold  a card  iaage. 

(Input)  Unit  no,  for  SPECS  f ile . 

(Ir.nut)  Urit  nc.  for  scratch  file. 

(Input)  Coded  version  and  date  cf  MIhOS. 
(Output)  will  be  negative  if  an  EOF  occurred 
on  file  ISPECS  before  a "PEGIN"  card  was  found. 
Otherwise,  returns  the  no.  of  card  images  that 
certained  useful  infornaticn. 


1 


I 


.1 


flethod : 

Since  the  scratch  file  is  also  used  for  row  and  coluan  naaes,  we  are 
stuck  with  the  sam®  Logical  Rfcord  Length,  b characters.  Four  logical 
records  will  l:ie  output  to  the  scratch  file  for  each  neaninaful  SPFCS 
card.  These  will  be  the  fcllcwing: 


KEYWORD 

IDENTIFIER 

INTEGEB 

FEAl 


(3  characters  ♦ 5 blanks) 

(8A  1,  tc  be  re-read  as  2X4) 
)6A1,  tc  te  re-tead  18) 
(8»,  1,  tc  be  re-tead  ns  E8.1) 


Depending  on  the  KEYWOPH,  any  of  the  IDPNTIhJSB,  INTEGfR  or  REAL  eay 
be  S blanks. 


1,  Read  cards  frea  file  ISFFCS  until  oa«  is  found  that  contains  "EEC" 
as  ehe  first  non-blank  characters.  Exit  with  IPECS  * -1  if  Er.d  cf 
File  occurs  first. 


2.  (BEGIN  found)  opcode  IVEFSN  and  output  heaiinn. 


I 


i 


I 

I 


3.  Pcoc«ss  «aci<  card  as  fcllovs. 

(a)  Pead  and  list  tb-i  card  iaaqe. 

(b)  Scan  to  the  first  non-blark  cbaracter.  Skip  card  if  coaplately 
blank,  or  if  first  character  is  •••*’,  which  indicates  a ccrsent. 

(c)  tr.it  if  the  hEYWCFD  is  "EPC". 

(d)  Taka  as  fCETWORD  the  3 characters  starting  froa,  and  including, 
the  first  non-tlank  character. 

(e)  Scan  to  the  neit  tlank  or 

(f)  Scan  past  all  blanks  cr  "=''s. 

(q)  The  next  character  dcteriinas  whether  the  second  itea  on  the 
card  is  ar.  identifier  or  a nushar.  If  it  is  cne  of  a list  of 
characters  contained  in  the  local  array  LDIGIT,  then  extract  up 
to  8 such  characters  and  treat  either  as  an  INTfSEP,  or  (if  it 
contains  a or  "E"  or  "C")  as  a FEM.  otherwise,  extract 

exactly  R characters  •■o  be  output  as  IDENTIFIEP. 

(h)  If  an  INTEGER  or  FEAL  has  cot  yet  been  found,  return  to  (e) . 

(i)  Output  the  U itess  KEYNCBC,  ICENTIFIt'E,  INTEGER  and  PEAL  to  the 
scratch  file. 

4.  If  a SPECS  file  was  found,  close  and  rewind  file  I5CB. 


i 

I 

!;i 


i 


LI 


I 


!■; 


I 
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SOBROOTlllE  SPECS2(  N SP EC , HPC H S , MCOLS  , MF LWS , 
1 n,M2,!!R,KNK,K,NE,NF1,ltF,NN,NN0,*!fl:(3,t’;s  ) 


I 


Subroatines  calletl: 


Called  by: 


None 


niNOs 


Parpcsa: 

To  aake  a second  pass  at  decodinq  inforiation  froa  the  SPFCS  file. 
The  scratch  file  v.nv  contains  a sequence  of  PFYilORDS  and  VHUSS  in 
fixed  fcteat.  It  reaains  tc  detaiiBice  lihich  KEYWORDS  they  are. 


l; 

i; 

1 


Par iiaters: 

INTEGE6  NSPFC  (Input)  Ihe  no.  of  KEYWORDS  to  be  read  froB 

the  scratch  file. 

The  re«aining  parameters  are  output.  They  are  the  various  INTEGET 
variables  defining  the  diiensicne  cf  the  problea  about  to  be  solved. 


Bet  hod : 

1.  Set  lost  paraneters  tc  their  default  values.  For  the  others,  it  is 
iupor*-ctnt  to  Know  whether  the  us‘'r  attempted  tc  define  values,  or 
not.  This  is  (i^tferminel  by  setting  ♦he  parameters  to  illegal 
values,  e.g,  0 or  -1,  and  then  testing  for  such  values  alter  all 
KEYWORDS  have  been  piccsssed. 

2.  Per  each  KBTWOFD,  read  the  variables  KEY,  ID,  N'JM,  DNDP  using 
foraatf  A3/,  2A4/,  18/,  E8.1  ). 

3.  Determine  which  Keyword  is  contained  in  KEY  by  searching  a list  of 
legal  Keywords.  The  use  cf  ID,  NOB  and  DNDK  depends  on  KEY,  Tn 
nest  cases,  only  one  of  It,  NUB,  or  DNDB  is  used. 

u.  when  all  Keywords  have  been  processed,  assign  default  values  to 
those  that  ware  net  specified. 

Print  a suaaary  of  the  paraaetecs  about  to  be  used,  and  exit. 


i 

li 

I 


I 


J 


V 
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STIBFCUTIIIE  SQUFI7(  H , N , N F ,NK  , N NO  , HS  , NS  , 

1 HE, HS, BF, BU, F, X,XN,G,GFC,TCLX,PTNOPM, PONOFM, DELETE  ) 

Subroutin«s  called:  Called  by: 

TELCOL  BGITN 


purpose ; 

To  rp«o»e  any  supertasics  that  are  very  close  to  a hound  and  do  not 
have  <1  aiqnif lean*  reduced  gradient  driving  theia  away  ^roa  th=it  riourd. 
fWill  tfe  called  CDly  after  ccnsttelned  s<-eps  ii  BGITN.) 


Paraieters: 

lOLX  (Input)  The  aeasure  of  closeness  to  a bound. 

P£Al»B  IINOa*!  (Input)  Nor*  of  PI,  to  scal^^  tht  tolerance  on 

reduced-qradient  coapenents. 

LOCICJIL  DBlEtP  (Irout)  To  be  passed  to  CLLCOL. 


Kethod : 

Suoerbdsicb  are  spli*'  up  depending  on  whether  *:heir  reduced  gradient 
is  s*8ll  or  no*-.  Ceding  could  ta  siiplified  if  ♦^he  putpost^  .stated 
above  is  followed  directly. 

If  a variable  is  tc  be  sgucezed,  its  state,  HS  (K) , is  aade  nonbdsic  at 
the  aprropriate  bound.  Arrays  G and  GFD  are  compacted,  and  DSLCOL 
does  tbfc'  same  for  HP,  x and  P. 

Throughout,  .superbasics  ate  processed  in  reverse  order  to  reduce  werk 
slightly  if  aore  than  one  gets  sgueered. 
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StJPPOOTINE  STATE  ( K^ISTATE  ) 

Subccatines  called: 

Mono 


Called  by: 
tflIVPF 


Purpose: 

To  load  appropriate  character  strings  into  the  array  TSTATE,  for 
output  with  a saved  bit  nap  or  the  final  sclution. 


Paraieters: 

INTEGER  K 

INTEGER  ISTATE{3) 


(Input)  Specifies  the  current  state  of  the 
picklen. 

(Output)  Will  contain  one  of  the  following 
strings,  4 characters  per  word: 

K-0  PBOCF'CING 

1 OPIIKAL  SCLN 

2 INFBAEIELE 

3 HNPCUNDED 

4 EXCESt'  ITSS 

5 ERROR  CCSDN 


Vi 


Het  hod : 
Trivial. 


SUBROUTINE  TP NS  VI ( N , K 1 , N E , IF , IE N , IPI ?E H , PC, R V , PREO EC, AC F ,OUT ,NR E J) 


Sabrcutines  called;  Callad  by; 

None  F4 


Purpcse; 


Given  t be  row  numbers  of  the  ncr.zeros  in  each  coluan  of  ^ square 
sparse  aatrix  A,  t^.is  rou*^,ine  finds  a ccluan  pet  mutation  P £:uch  that 
AF  has  a nonzero  diaqonal.  If  this  is  rot  possitl*',  a jerautation  is 
found  such  ‘hat  ♦ho  no.  of  consecutive  r.crzero  diagonals  in  is  as 
larqe  as  possible. 


Paraaeters; 


INTEG-^F 

INTSGSF 

INTEGER 

IVTPGi:F*2 

INTEGER* 2 
INTEGEF*2 


INTKGPF*2 
INr£GEh*2 
INTEGSR*2 
INTRGFR*2 
INTEGER *2 
INTEGER 


N 

PI 

NZ 

IP  (N1) 

IHN  fNZ) 
IPIVRW  (N) 


PC(N) 

PV  (N) 
PEEOHD (N) 
ACP  (N) 

CUT  (N< 
NhEJ 


ftnput ) 
(Input) 
(Input) 
(Input) 
colun  D . 
(Input) 


Dimansicn  of  the  matrix. 

N*1. 

No.  cf  noDzeros  in  the  matrix, 
Points  to  the  firs*  nc:.zeto 


in  each 


List  of  row  numbers  for  nonzero?. 


(Output)  Will  contain  the  requirt d permutation. 


as  an  N *N  ar i a y of 
elements  A (I,  I r I VFW  (i) ) 
0 ; see 


Ccnsidering  the  matrix 
elements  A(I,J),  the 
will  be  nonzero  (except  if  IPIVPW  (I) 

NHEJ  below) . 
fcork  vector. 

Work  vector. 

Work  vector, 
work  vector. 

Work  vector. 

(Output)  will  be  nonzero  if  the  matrix  is 
structurallv  slnqular  (i.e.  singalar  for  all 
pcFsible  values  of  the  nonzeros).  Certain 
elener*?  of  IPIVPW  vil'  then  be  zero  and  it  is 
possible  *(*^CQfore  tc  detornine  which  coluan.s 
cf  the  matrix  most  be  replaced  to  obtain  a 
nonsirqulac  matrix. 


nathod : 

This  routine  Is  derived  directly  froi  Harwell  subroutine  (!C'?6B  by  I.S. 
Duff  (dated  Doc  19V h)  , which  cay  rirce  hive  replaced  by  NCIIA. 


SnBFOPTINE  OVFACK(  fl , NI , MP 1 , H A , H E , A , Y , JQ  ) 


SubcoQtines  called;  Called  by; 

None  ADDCCL  DPIVER  FACTOR  LPITN 

Fr;iT.J  SOLM 


Purpose ; 

Tc  expand  the  packed  cclusn  A (JO)  into  the  vector  Y. 


Net hod : 

The  paranieters  and  »ethod  are  cbvious  from  the  code.  Kcte 
negat3.ve  tew  idices  HA  (J)  nay  occur  in  special  applications 
(lINCS/GRf; ) , but  only  their  atsclote  value  is  relevant  hare. 


SUBFOOTINE  UNPACK  ( N , »E , N P 1 , H A , H E , A , Y , JQ  ) 
IN-^EGEB*?  VA  (NE)  ,IF  (NF1) 

PEAL  A(NE) 

PEA1*8  Y(?H 

C 

C UNPACK  THE  JQ-1H  CCLUHN  OF  A 

C 

FO  10  I»1,!! 

10  Y (I)  « 0.0 
C 

J1  » HZ(JO)  +1 
J2  * HE(JQ*1) 

C 

DC  20  1*J1,J2 
IB  * HA  (I) 

IB  ’ IADS  (IR) 

T(IH)  = A{H 
20  CCNTINUE 
C 

BETURN 

C END  OF  UNPACK 

END 
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11.  PERFORMANCE 


The  proving  ground  for  MINOS  has  been  the  U.S.  energy-economic  mooels 
of  Alan  Manne.  These  are  the  models  ETA  (Manne,  1976)  and  ETA-MACRO 
(Manne,  1977),  of  which  several  hundred  cases  have  been  run  during  1975-77. 

In  ETA-MACRO,  the  constraint  matrix  is  approximately  430  x 750  with 
2800  nonzero  coefficients.  The  objective  function  involves  80  variables 
nonli nearly.  It  is  of  the  form 

f(K,E,N,I,EC)  = I 6 log  C 
t=i 

whprp 

C,  > y,  - (ac^  + - EC^  . 

= (K.  - , 

d.  = (E,  . 

and  the  unknowns  are  the  five  sets  of  variables  E^,  N^,  I^,  EC^. 

Clearly  log  C^  is  a non-separable  function  of  these  variables. 

On  an  IBM  System  370  the  problem  can  be  run  i*^  a 384K  region.  This 
core  requirement  is  made  up  approximately  as  follows: 

Program  code  184K 

REAL*8  7(20000)  160K 

I/O  buffers  40K 

384 K bytes. 

From  a cold  start,  using  CRASH  OPTION  1 with  no  special  initialization 
of  the  nonlinear  variables,  the  run  time  on  an  IBM  370/168  is  typically 
about  100  seconds  for  1100  iterations.  (The  first  500  iterations  are  primal 
simplex  iterations  obtaining  an  initial  feasible  solution.)  Restarts  with 
different  bound  sets  require  a few  hundred  further  iterations. 

Several  factors  contribute  to  the  efficiency  of  these  runs.  The  con- 
straints are  very  sparse,  and  the  LU  factorization  of  each  basis  typically 
has  4 bumps  and  only  8 spikes.  The  number  of  superbasic  variables,  and 
hence  the  dimension  of  the  reduced  Hessian  approximation,  remains  at  around 
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30,  which  is  well  below  the  theoretical  upper  bound  of  80.  Finally, 
although  a relatively  accurate  search  is  used  (ETASCH  = 0 01),  the  line- 


search  procedure  of  Gill  et  al  (1976)  in  this 
and  gradient  evaluations  per  iteration. 


case  averages  only  2 function 


12.  ADDENDA  TO  MINOS  USER'S  GUIDE 


The  User's  Guide  accurately  describes  MINOS  Version  3.1.  Some  of  the 
changes  incorporated  in  MINOS  Version  3.3  will  be  apparent  to  the  user,  as 
follows: 

1.  Two  additional  parameters  may  be  specified  in  the  SPECS  file,  by  cards 
of  the  form 

DO  TOLERANCE  l.OE-5 

FEASIBILITY  TCL  l.OE-4 

These  reset  the  variables  TOLOJS  and  TOLX  respectively  (see  section  5). 

2.  The  default  value  for  the  feasibility  tolerance  TOLX  is  now 
max(v^,  10'^). 

3.  Subroutines  ALIGN  and  SQUEEZ  have  been  deleted.  The  keyword 
ALIGNMENT  TOLFRANCE  no  longer  has  the  effect  described  on  p.  40  of  the 
User's  Guide. 

4.  Variables  that  are  specified  to  be  stiperbasic  in 

• an  INITIAL  bounds  set  in  the  MPS  file 

• an  OLD  BASIS  FILE 

will  be  initialized  at  the  specified  values  regardless  of  feasibility. 
(Exceptions:  "infinite"  values  are  replaced  by  zero;  if  the  SUPERBASICS 
LIMIT  has  been  reached,  variables  will  be  made  nonbasic  at  the  nearest 
bound.)  This  facilitates  initialization  at  intentionally  infeasible 
points,  and  exact  reconstruction  of  solutions  from  previous  runs.  The 
following  parts  of  the  User's  Guide  should  be  altered  accordingly: 
p.  29,  8(a). 

p.  62,  8. 

5.  MPS  FILE  is  a synonym  for  INPUT  FILE. 
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